Differential D5594 Diff 18104 extern/quadriflow/3rd/lemon-1.3.1/lemon/grosso_locatelli_pullan_mc.h

Changeset View

Standalone View

extern/quadriflow/3rd/lemon-1.3.1/lemon/grosso_locatelli_pullan_mc.h

This file was added.

				/* -- mode: C++; indent-tabs-mode: nil; --
				*
				* This file is a part of LEMON, a generic C++ optimization library.
				*
				* Copyright (C) 2003-2013
				* Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
				* (Egervary Research Group on Combinatorial Optimization, EGRES).
				*
				* Permission to use, modify and distribute this software is granted
				* provided that this copyright notice appears in all copies. For
				* precise terms see the accompanying LICENSE file.
				*
				* This software is provided "AS IS" with no warranty of any kind,
				* express or implied, and with no claim as to its suitability for any
				* purpose.
				*
				*/

				#ifndef LEMON_GROSSO_LOCATELLI_PULLAN_MC_H
				#define LEMON_GROSSO_LOCATELLI_PULLAN_MC_H

				/// \ingroup approx_algs
				///
				/// \file
				/// \brief The iterated local search algorithm of Grosso, Locatelli, and Pullan
				/// for the maximum clique problem

				#include <vector>
				#include <limits>
				#include <lemon/core.h>
				#include <lemon/random.h>

				namespace lemon {

				/// \addtogroup approx_algs
				/// @{

				/// \brief Implementation of the iterated local search algorithm of Grosso,
				/// Locatelli, and Pullan for the maximum clique problem
				///
				/// \ref GrossoLocatelliPullanMc implements the iterated local search
				/// algorithm of Grosso, Locatelli, and Pullan for solving the \e maximum
				/// \e clique \e problem \cite grosso08maxclique.
				/// It is to find the largest complete subgraph (\e clique) in an
				/// undirected graph, i.e., the largest set of nodes where each
				/// pair of nodes is connected.
				///
				/// This class provides a simple but highly efficient and robust heuristic
				/// method that quickly finds a quite large clique, but not necessarily the
				/// largest one.
				/// The algorithm performs a certain number of iterations to find several
				/// cliques and selects the largest one among them. Various limits can be
				/// specified to control the running time and the effectiveness of the
				/// search process.
				///
				/// \tparam GR The undirected graph type the algorithm runs on.
				///
				/// \note %GrossoLocatelliPullanMc provides three different node selection
				/// rules, from which the most powerful one is used by default.
				/// For more information, see \ref SelectionRule.
				template <typename GR>
				class GrossoLocatelliPullanMc
				{
				public:

				/// \brief Constants for specifying the node selection rule.
				///
				/// Enum type containing constants for specifying the node selection rule
				/// for the \ref run() function.
				///
				/// During the algorithm, nodes are selected for addition to the current
				/// clique according to the applied rule.
				/// In general, the PENALTY_BASED rule turned out to be the most powerful
				/// and the most robust, thus it is the default option.
				/// However, another selection rule can be specified using the \ref run()
				/// function with the proper parameter.
				enum SelectionRule {

				/// A node is selected randomly without any evaluation at each step.
				RANDOM,

				/// A node of maximum degree is selected randomly at each step.
				DEGREE_BASED,

				/// A node of minimum penalty is selected randomly at each step.
				/// The node penalties are updated adaptively after each stage of the
				/// search process.
				PENALTY_BASED
				};

				/// \brief Constants for the causes of search termination.
				///
				/// Enum type containing constants for the different causes of search
				/// termination. The \ref run() function returns one of these values.
				enum TerminationCause {

				/// The iteration count limit is reached.
				ITERATION_LIMIT,

				/// The step count limit is reached.
				STEP_LIMIT,

				/// The clique size limit is reached.
				SIZE_LIMIT
				};

				private:

				TEMPLATE_GRAPH_TYPEDEFS(GR);

				typedef std::vector<int> IntVector;
				typedef std::vector<char> BoolVector;
				typedef std::vector<BoolVector> BoolMatrix;
				// Note: vector<char> is used instead of vector<bool> for efficiency reasons

				// The underlying graph
				const GR &_graph;
				IntNodeMap _id;

				// Internal matrix representation of the graph
				BoolMatrix _gr;
				int _n;

				// Search options
				bool _delta_based_restart;
				int _restart_delta_limit;

				// Search limits
				int _iteration_limit;
				int _step_limit;
				int _size_limit;

				// The current clique
				BoolVector _clique;
				int _size;

				// The best clique found so far
				BoolVector _best_clique;
				int _best_size;

				// The "distances" of the nodes from the current clique.
				// _delta[u] is the number of nodes in the clique that are
				// not connected with u.
				IntVector _delta;

				// The current tabu set
				BoolVector _tabu;

				// Random number generator
				Random _rnd;

				private:

				// Implementation of the RANDOM node selection rule.
				class RandomSelectionRule
				{
				private:

				// References to the algorithm instance
				const BoolVector &_clique;
				const IntVector &_delta;
				const BoolVector &_tabu;
				Random &_rnd;

				// Pivot rule data
				int _n;

				public:

				// Constructor
				RandomSelectionRule(GrossoLocatelliPullanMc &mc) :
				_clique(mc._clique), _delta(mc._delta), _tabu(mc._tabu),
				_rnd(mc._rnd), _n(mc._n)
				{}

				// Return a node index for a feasible add move or -1 if no one exists
				int nextFeasibleAddNode() const {
				int start_node = _rnd[_n];
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0 && !_tabu[i]) return i;
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0 && !_tabu[i]) return i;
				}
				return -1;
				}

				// Return a node index for a feasible swap move or -1 if no one exists
				int nextFeasibleSwapNode() const {
				int start_node = _rnd[_n];
				for (int i = start_node; i != _n; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i]) return i;
				}
				for (int i = 0; i != start_node; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i]) return i;
				}
				return -1;
				}

				// Return a node index for an add move or -1 if no one exists
				int nextAddNode() const {
				int start_node = _rnd[_n];
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0) return i;
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0) return i;
				}
				return -1;
				}

				// Update internal data structures between stages (if necessary)
				void update() {}

				}; //class RandomSelectionRule


				// Implementation of the DEGREE_BASED node selection rule.
				class DegreeBasedSelectionRule
				{
				private:

				// References to the algorithm instance
				const BoolVector &_clique;
				const IntVector &_delta;
				const BoolVector &_tabu;
				Random &_rnd;

				// Pivot rule data
				int _n;
				IntVector _deg;

				public:

				// Constructor
				DegreeBasedSelectionRule(GrossoLocatelliPullanMc &mc) :
				_clique(mc._clique), _delta(mc._delta), _tabu(mc._tabu),
				_rnd(mc._rnd), _n(mc._n), _deg(_n)
				{
				for (int i = 0; i != _n; i++) {
				int d = 0;
				BoolVector &row = mc._gr[i];
				for (int j = 0; j != _n; j++) {
				if (row[j]) d++;
				}
				_deg[i] = d;
				}
				}

				// Return a node index for a feasible add move or -1 if no one exists
				int nextFeasibleAddNode() const {
				int start_node = _rnd[_n];
				int node = -1, max_deg = -1;
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0 && !_tabu[i] && _deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0 && !_tabu[i] && _deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				return node;
				}

				// Return a node index for a feasible swap move or -1 if no one exists
				int nextFeasibleSwapNode() const {
				int start_node = _rnd[_n];
				int node = -1, max_deg = -1;
				for (int i = start_node; i != _n; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i] &&
				_deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i] &&
				_deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				return node;
				}

				// Return a node index for an add move or -1 if no one exists
				int nextAddNode() const {
				int start_node = _rnd[_n];
				int node = -1, max_deg = -1;
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0 && _deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0 && _deg[i] > max_deg) {
				node = i;
				max_deg = _deg[i];
				}
				}
				return node;
				}

				// Update internal data structures between stages (if necessary)
				void update() {}

				}; //class DegreeBasedSelectionRule


				// Implementation of the PENALTY_BASED node selection rule.
				class PenaltyBasedSelectionRule
				{
				private:

				// References to the algorithm instance
				const BoolVector &_clique;
				const IntVector &_delta;
				const BoolVector &_tabu;
				Random &_rnd;

				// Pivot rule data
				int _n;
				IntVector _penalty;

				public:

				// Constructor
				PenaltyBasedSelectionRule(GrossoLocatelliPullanMc &mc) :
				_clique(mc._clique), _delta(mc._delta), _tabu(mc._tabu),
				_rnd(mc._rnd), _n(mc._n), _penalty(_n, 0)
				{}

				// Return a node index for a feasible add move or -1 if no one exists
				int nextFeasibleAddNode() const {
				int start_node = _rnd[_n];
				int node = -1, min_p = std::numeric_limits<int>::max();
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0 && !_tabu[i] && _penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0 && !_tabu[i] && _penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				return node;
				}

				// Return a node index for a feasible swap move or -1 if no one exists
				int nextFeasibleSwapNode() const {
				int start_node = _rnd[_n];
				int node = -1, min_p = std::numeric_limits<int>::max();
				for (int i = start_node; i != _n; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i] &&
				_penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (!_clique[i] && _delta[i] == 1 && !_tabu[i] &&
				_penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				return node;
				}

				// Return a node index for an add move or -1 if no one exists
				int nextAddNode() const {
				int start_node = _rnd[_n];
				int node = -1, min_p = std::numeric_limits<int>::max();
				for (int i = start_node; i != _n; i++) {
				if (_delta[i] == 0 && _penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				for (int i = 0; i != start_node; i++) {
				if (_delta[i] == 0 && _penalty[i] < min_p) {
				node = i;
				min_p = _penalty[i];
				}
				}
				return node;
				}

				// Update internal data structures between stages (if necessary)
				void update() {}

				}; //class PenaltyBasedSelectionRule

				public:

				/// \brief Constructor.
				///
				/// Constructor.
				/// The global \ref rnd "random number generator instance" is used
				/// during the algorithm.
				///
				/// \param graph The undirected graph the algorithm runs on.
				GrossoLocatelliPullanMc(const GR& graph) :
				_graph(graph), _id(_graph), _rnd(rnd)
				{
				initOptions();
				}

				/// \brief Constructor with random seed.
				///
				/// Constructor with random seed.
				///
				/// \param graph The undirected graph the algorithm runs on.
				/// \param seed Seed value for the internal random number generator
				/// that is used during the algorithm.
				GrossoLocatelliPullanMc(const GR& graph, int seed) :
				_graph(graph), _id(_graph), _rnd(seed)
				{
				initOptions();
				}

				/// \brief Constructor with random number generator.
				///
				/// Constructor with random number generator.
				///
				/// \param graph The undirected graph the algorithm runs on.
				/// \param random A random number generator that is used during the
				/// algorithm.
				GrossoLocatelliPullanMc(const GR& graph, const Random& random) :
				_graph(graph), _id(_graph), _rnd(random)
				{
				initOptions();
				}

				/// \name Execution Control
				/// The \ref run() function can be used to execute the algorithm.\n
				/// The functions \ref iterationLimit(int), \ref stepLimit(int), and
				/// \ref sizeLimit(int) can be used to specify various limits for the
				/// search process.

				/// @{

				/// \brief Sets the maximum number of iterations.
				///
				/// This function sets the maximum number of iterations.
				/// Each iteration of the algorithm finds a maximal clique (but not
				/// necessarily the largest one) by performing several search steps
				/// (node selections).
				///
				/// This limit controls the running time and the success of the
				/// algorithm. For larger values, the algorithm runs slower, but it more
				/// likely finds larger cliques. For smaller values, the algorithm is
				/// faster but probably gives worse results.
				///
				/// The default value is \c 1000.
				/// \c -1 means that number of iterations is not limited.
				///
				/// \warning You should specify a reasonable limit for the number of
				/// iterations and/or the number of search steps.
				///
				/// \return <tt>(*this)</tt>
				///
				/// \sa stepLimit(int)
				/// \sa sizeLimit(int)
				GrossoLocatelliPullanMc& iterationLimit(int limit) {
				_iteration_limit = limit;
				return *this;
				}

				/// \brief Sets the maximum number of search steps.
				///
				/// This function sets the maximum number of elementary search steps.
				/// Each iteration of the algorithm finds a maximal clique (but not
				/// necessarily the largest one) by performing several search steps
				/// (node selections).
				///
				/// This limit controls the running time and the success of the
				/// algorithm. For larger values, the algorithm runs slower, but it more
				/// likely finds larger cliques. For smaller values, the algorithm is
				/// faster but probably gives worse results.
				///
				/// The default value is \c -1, which means that number of steps
				/// is not limited explicitly. However, the number of iterations is
				/// limited and each iteration performs a finite number of search steps.
				///
				/// \warning You should specify a reasonable limit for the number of
				/// iterations and/or the number of search steps.
				///
				/// \return <tt>(*this)</tt>
				///
				/// \sa iterationLimit(int)
				/// \sa sizeLimit(int)
				GrossoLocatelliPullanMc& stepLimit(int limit) {
				_step_limit = limit;
				return *this;
				}

				/// \brief Sets the desired clique size.
				///
				/// This function sets the desired clique size that serves as a search
				/// limit. If a clique of this size (or a larger one) is found, then the
				/// algorithm terminates.
				///
				/// This function is especially useful if you know an exact upper bound
				/// for the size of the cliques in the graph or if any clique above
				/// a certain size limit is sufficient for your application.
				///
				/// The default value is \c -1, which means that the size limit is set to
				/// the number of nodes in the graph.
				///
				/// \return <tt>(*this)</tt>
				///
				/// \sa iterationLimit(int)
				/// \sa stepLimit(int)
				GrossoLocatelliPullanMc& sizeLimit(int limit) {
				_size_limit = limit;
				return *this;
				}

				/// \brief The maximum number of iterations.
				///
				/// This function gives back the maximum number of iterations.
				/// \c -1 means that no limit is specified.
				///
				/// \sa iterationLimit(int)
				int iterationLimit() const {
				return _iteration_limit;
				}

				/// \brief The maximum number of search steps.
				///
				/// This function gives back the maximum number of search steps.
				/// \c -1 means that no limit is specified.
				///
				/// \sa stepLimit(int)
				int stepLimit() const {
				return _step_limit;
				}

				/// \brief The desired clique size.
				///
				/// This function gives back the desired clique size that serves as a
				/// search limit. \c -1 means that this limit is set to the number of
				/// nodes in the graph.
				///
				/// \sa sizeLimit(int)
				int sizeLimit() const {
				return _size_limit;
				}

				/// \brief Runs the algorithm.
				///
				/// This function runs the algorithm. If one of the specified limits
				/// is reached, the search process terminates.
				///
				/// \param rule The node selection rule. For more information, see
				/// \ref SelectionRule.
				///
				/// \return The termination cause of the search. For more information,
				/// see \ref TerminationCause.
				TerminationCause run(SelectionRule rule = PENALTY_BASED)
				{
				init();
				switch (rule) {
				case RANDOM:
				return start<RandomSelectionRule>();
				case DEGREE_BASED:
				return start<DegreeBasedSelectionRule>();
				default:
				return start<PenaltyBasedSelectionRule>();
				}
				}

				/// @}

				/// \name Query Functions
				/// The results of the algorithm can be obtained using these functions.\n
				/// The run() function must be called before using them.

				/// @{

				/// \brief The size of the found clique
				///
				/// This function returns the size of the found clique.
				///
				/// \pre run() must be called before using this function.
				int cliqueSize() const {
				return _best_size;
				}

				/// \brief Gives back the found clique in a \c bool node map
				///
				/// This function gives back the characteristic vector of the found
				/// clique in the given node map.
				/// It must be a \ref concepts::WriteMap "writable" node map with
				/// \c bool (or convertible) value type.
				///
				/// \pre run() must be called before using this function.
				template <typename CliqueMap>
				void cliqueMap(CliqueMap &map) const {
				for (NodeIt n(_graph); n != INVALID; ++n) {
				map[n] = static_cast<bool>(_best_clique[_id[n]]);
				}
				}

				/// \brief Iterator to list the nodes of the found clique
				///
				/// This iterator class lists the nodes of the found clique.
				/// Before using it, you must allocate a GrossoLocatelliPullanMc instance
				/// and call its \ref GrossoLocatelliPullanMc::run() "run()" method.
				///
				/// The following example prints out the IDs of the nodes in the found
				/// clique.
				/// \code
				/// GrossoLocatelliPullanMc<Graph> mc(g);
				/// mc.run();
				/// for (GrossoLocatelliPullanMc<Graph>::CliqueNodeIt n(mc);
				/// n != INVALID; ++n)
				/// {
				/// std::cout << g.id(n) << std::endl;
				/// }
				/// \endcode
				class CliqueNodeIt
				{
				private:
				NodeIt _it;
				BoolNodeMap _map;

				public:

				/// Constructor

				/// Constructor.
				/// \param mc The algorithm instance.
				CliqueNodeIt(const GrossoLocatelliPullanMc &mc)
				: _map(mc._graph)
				{
				mc.cliqueMap(_map);
				for (_it = NodeIt(mc._graph); _it != INVALID && !_map[_it]; ++_it) ;
				}

				/// Conversion to \c Node
				operator Node() const { return _it; }

				bool operator==(Invalid) const { return _it == INVALID; }
				bool operator!=(Invalid) const { return _it != INVALID; }

				/// Next node
				CliqueNodeIt &operator++() {
				for (++_it; _it != INVALID && !_map[_it]; ++_it) ;
				return *this;
				}

				/// Postfix incrementation

				/// Postfix incrementation.
				///
				/// \warning This incrementation returns a \c Node, not a
				/// \c CliqueNodeIt as one may expect.
				typename GR::Node operator++(int) {
				Node n=*this;
				++(*this);
				return n;
				}

				};

				/// @}

				private:

				// Initialize search options and limits
				void initOptions() {
				// Search options
				_delta_based_restart = true;
				_restart_delta_limit = 4;

				// Search limits
				_iteration_limit = 1000;
				_step_limit = -1; // this is disabled by default
				_size_limit = -1; // this is disabled by default
				}

				// Adds a node to the current clique
				void addCliqueNode(int u) {
				if (_clique[u]) return;
				_clique[u] = true;
				_size++;
				BoolVector &row = _gr[u];
				for (int i = 0; i != _n; i++) {
				if (!row[i]) _delta[i]++;
				}
				}

				// Removes a node from the current clique
				void delCliqueNode(int u) {
				if (!_clique[u]) return;
				_clique[u] = false;
				_size--;
				BoolVector &row = _gr[u];
				for (int i = 0; i != _n; i++) {
				if (!row[i]) _delta[i]--;
				}
				}

				// Initialize data structures
				void init() {
				_n = countNodes(_graph);
				int ui = 0;
				for (NodeIt u(_graph); u != INVALID; ++u) {
				_id[u] = ui++;
				}
				_gr.clear();
				_gr.resize(_n, BoolVector(_n, false));
				ui = 0;
				for (NodeIt u(_graph); u != INVALID; ++u) {
				for (IncEdgeIt e(_graph, u); e != INVALID; ++e) {
				int vi = _id[_graph.runningNode(e)];
				_gr[ui][vi] = true;
				_gr[vi][ui] = true;
				}
				++ui;
				}

				_clique.clear();
				_clique.resize(_n, false);
				_size = 0;
				_best_clique.clear();
				_best_clique.resize(_n, false);
				_best_size = 0;
				_delta.clear();
				_delta.resize(_n, 0);
				_tabu.clear();
				_tabu.resize(_n, false);
				}

				// Executes the algorithm
				template <typename SelectionRuleImpl>
				TerminationCause start() {
				if (_n == 0) return SIZE_LIMIT;
				if (_n == 1) {
				_best_clique[0] = true;
				_best_size = 1;
				return SIZE_LIMIT;
				}

				// Iterated local search algorithm
				const int max_size = _size_limit >= 0 ? _size_limit : _n;
				const int max_restart = _iteration_limit >= 0 ?
				_iteration_limit : std::numeric_limits<int>::max();
				const int max_select = _step_limit >= 0 ?
				_step_limit : std::numeric_limits<int>::max();

				SelectionRuleImpl sel_method(*this);
				int select = 0, restart = 0;
				IntVector restart_nodes;
				while (select < max_select && restart < max_restart) {

				// Perturbation/restart
				restart++;
				if (_delta_based_restart) {
				restart_nodes.clear();
				for (int i = 0; i != _n; i++) {
				if (_delta[i] >= _restart_delta_limit)
				restart_nodes.push_back(i);
				}
				}
				int rs_node = -1;
				if (restart_nodes.size() > 0) {
				rs_node = restart_nodes[_rnd[restart_nodes.size()]];
				} else {
				rs_node = _rnd[_n];
				}
				BoolVector &row = _gr[rs_node];
				for (int i = 0; i != _n; i++) {
				if (_clique[i] && !row[i]) delCliqueNode(i);
				}
				addCliqueNode(rs_node);

				// Local search
				_tabu.clear();
				_tabu.resize(_n, false);
				bool tabu_empty = true;
				int max_swap = _size;
				while (select < max_select) {
				select++;
				int u;
				if ((u = sel_method.nextFeasibleAddNode()) != -1) {
				// Feasible add move
				addCliqueNode(u);
				if (tabu_empty) max_swap = _size;
				}
				else if ((u = sel_method.nextFeasibleSwapNode()) != -1) {
				// Feasible swap move
				int v = -1;
				BoolVector &row = _gr[u];
				for (int i = 0; i != _n; i++) {
				if (_clique[i] && !row[i]) {
				v = i;
				break;
				}
				}
				addCliqueNode(u);
				delCliqueNode(v);
				_tabu[v] = true;
				tabu_empty = false;
				if (--max_swap <= 0) break;
				}
				else if ((u = sel_method.nextAddNode()) != -1) {
				// Non-feasible add move
				addCliqueNode(u);
				}
				else break;
				}
				if (_size > _best_size) {
				_best_clique = _clique;
				_best_size = _size;
				if (_best_size >= max_size) return SIZE_LIMIT;
				}
				sel_method.update();
				}

				return (restart >= max_restart ? ITERATION_LIMIT : STEP_LIMIT);
				}

				}; //class GrossoLocatelliPullanMc

				///@}

				} //namespace lemon

				#endif //LEMON_GROSSO_LOCATELLI_PULLAN_MC_H