/* Copyright (c) 1997-2006
   Ewgenij Gawrilow, Michael Joswig (Technische Universitaet Berlin, Germany)
   http://www.math.tu-berlin.de/polymake,  mailto:polymake@math.tu-berlin.de

   This program is free software; you can redistribute it and/or modify it
   under the terms of the GNU General Public License as published by the
   Free Software Foundation; either version 2, or (at your option) any
   later version: http://www.gnu.org/licenses/gpl.txt.

   This program is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.
*/

#ifndef _POLYMAKE_SHORTEST_PATH_BF_H
#define _POLYMAKE_SHORTEST_PATH_BF_H "$Project: polymake $$Id: shortest_path_bf.h 6920 2006-02-01 17:39:17Z gawrilow $"

#include <shortest_path_common.h>
#include <list>

namespace polymake { namespace graph {

template <typename _Adapter>
class shortest_path_bf
   : public shortest_path< shortest_path_bf<_Adapter>, _Adapter > {
protected:
   typedef shortest_path<shortest_path_bf> super;
   friend class shortest_path<shortest_path_bf>;

public:
   typedef typename _Adapter::node_id node_id;
   typedef typename _Adapter::edge_id edge_id;
   typedef typename _Adapter::distance_type distance_type;

   shortest_path_bf(const typename _Adapter::basic_graph_type& G)
      : super(G) { }

   shortest_path_bf(const _Adapter& G)
      : super(G) { }

protected:
   // first: the first child
   // second: the next sibling
   struct spanning_tree_node : std::pair<node_id, node_id> {
      spanning_tree_node()
	 : std::pair<node_id,node_id>(_Adapter::void_node_id(), _Adapter::void_node_id()) { }
   };
   typedef typename _Adapter::template property_map<node_id, spanning_tree_node>::type spanning_tree_type;
   typedef typename super::terminal_map visited_map_type;

   void init_sources(const node_id source)
   {
      super::init_sources(source);
      queue.push_back(source);
      visited[source]=true;
#ifdef POLY_DEBUG
      dump_source=source;
#endif
   }
   using shortest_path<shortest_path_bf,_Adapter>::init_sources;

   template <typename Targets>
   void init_targets(const Targets&) { }

   void init()
   {
      super::init();
      this->data.init(spanning_tree);
      this->data.init(visited);
      this->data.init(excess);
      queue.clear();
   }

   node_id nearest_target(const node_id target)
   {
      if (!visited[target])
	 throw std::runtime_error("shortest_path_bf: target unreachable from the source");
      return target;
   }

   template <typename Targets>
   node_id nearest_target(const Targets& targets)
   {
      node_id n=this->data.void_node_id();
      for (typename Entire<Targets>::const_iterator t=entire(targets);  !t.at_end();  ++t) {
	 if (visited[*t])
	    if (n==this->data.void_node_id() || this->dist[*t] < this->dist[n])
	       n=*t;
      }
      if (n==this->data.void_node_id())
	 throw std::runtime_error("shortest_path_bf: target unreachable from the source");
      return n;
   }

   node_id spanning_tree_next(node_id cur, node_id start)
   {
      node_id child=spanning_tree[cur].first;
      if (child != this->data.void_node_id())
	 return child;

      while (cur!=start) {
	 node_id sibling=spanning_tree[cur].second;
	 if (sibling != this->data.void_node_id())
	    return sibling;
	 cur=this->data.from_node(this->arc2s[cur]);
      }

      return start;
   }

#ifdef POLY_DEBUG
   node_id dump_source;

   void dump_subtree(node_id n)
   {
      cout << n;
      if (spanning_tree[n].first != data.void_node_id()) {
	 cout << "( ";
	 dump_subtree(spanning_tree[n].first);
	 cout << " )";
      }
      if (spanning_tree[n].second != data.void_node_id()) {
	 cout << ", ";
	 dump_subtree(spanning_tree[n].second);
      }
   }
   void dump_spanning_tree()
   {
      dump_subtree(dump_source);
   }
#endif

   void set_excess(distance_type& exc, const distance_type delta)
   {
      if (exc<0 || exc>delta) exc=delta;
   }

   bool reset_excess(distance_type& exc, const distance_type delta)
   {
      return exc>=0 && (exc+=delta)<0;
   }

   template <typename Targets>
   int do_find (const Targets& targets)
   {
      while (!queue.empty()) {
	 node_id n=queue.front();  queue.pop_front();

	 const distance_type d=this->dist[n];
	 excess[n]=distance_type(-1);
#ifdef POLY_DEBUG
	 cout << "visiting " << n << endl;
#endif
	 typename pm::attrib<typename _Adapter::OutEdges>::plus_const_ref arcs=this->data.out_edges(n);
	 for (typename Entire<typename _Adapter::OutEdges>::const_iterator ai=entire(arcs); !ai.at_end(); ++ai) {
	    edge_id a=*ai;
	    node_id t=this->data.to_node(a);
	    const distance_type dnew=d + this->data.get_edge_weight(a);

	    if (!visited[t]) {
	       queue.push_back(t);
	       visited[t]=true;
	       this->dist[t]=dnew;
	       this->arc2s[t]=a;
	       spanning_tree[t].second=spanning_tree[n].first;
	       spanning_tree[n].first=t;
#ifdef POLY_DEBUG
	       cout << "new " << t << "\ndist: " << dist.get_container() << "\n";
	       dump_spanning_tree();
	       cout << endl;
#endif
	    } else {
	       const distance_type dist_delta=dnew-this->dist[t];
	       if (dist_delta<0) {
		  this->dist[t]=dnew;
		  node_id old_parent=this->data.from_node(this->arc2s[t]), sibling=spanning_tree[old_parent].first, s;
		  this->arc2s[t]=a;
		  if (sibling==t) {
		     spanning_tree[old_parent].first=spanning_tree[t].second;
		  } else {
		     node_id s;
		     while ((s=spanning_tree[sibling].second) != t) sibling=s;
		     spanning_tree[sibling].second=spanning_tree[t].second;
		  }
		  spanning_tree[t].second=spanning_tree[n].first;
		  spanning_tree[n].first=t;

		  s=t;
		  while (true) {
		     if (reset_excess(excess[s], dist_delta))
			queue.push_front(s);
		     if ((s=spanning_tree_next(s,t)) == t)
			break;
		     if (s==n)
			throw std::runtime_error("shortest_path_bf: encountered negative cycle");
		     this->dist[s]+=dist_delta;
		  }

#ifdef POLY_DEBUG
		  cout << "updated " << t << "\ndist: " << dist.get_container() << "\n";
		  dump_spanning_tree();
		  cout << endl;
#endif
	       } else {
		  set_excess(excess[n], dist_delta);
	       }
	    }
	 }
#ifdef POLY_DEBUG
	 cout << "excess=" << excess[n] << endl;
#endif
      }
      return nearest_target(targets);
   }

   spanning_tree_type spanning_tree;
   visited_map_type visited;
   typename super::distance_map excess;
   std::list<node_id> queue;
};

} } // end namespace polymake::graph

#endif // _POLYMAKE_SHORTEST_PATH_BF_H

// Local Variables:
// mode:C++
// c-basic-offset:3
// End: