/*
 *  Main authors:
 *     Patrick Pekczynski <pekczynski@ps.uni-sb.de>
 *
 *  Copyright:
 *     Patrick Pekczynski, 2004
 *
 *  Last modified:
 *     $Date: 2006-08-04 16:03:26 +0200 (Fri, 04 Aug 2006) $ by $Author: schulte $
 *     $Revision: 3512 $
 *
 *  This file is part of Gecode, the generic constraint
 *  development environment:
 *     http://www.gecode.org
 *
 *  See the file "LICENSE" for information on usage and
 *  redistribution of this file, and for a
 *     DISCLAIMER OF ALL WARRANTIES.
 *
 */	

#include "gecode/support/static-stack.hh"

namespace Gecode { namespace Int { namespace Sortedness {

  /**
   *   \brief Compute the sccs of the oriented intersection-graph
   *
   *   An y-node \f$y_j\f$ and its corresponding matching mate
   *   \f$x_{\phi(j)}\f$ form the smallest possible scc, since both
   *   edges \f$e_1(y_j, x_{\phi(j)})\f$ and \f$e_2(x_{\phi(j)},y_j)\f$
   *   are both contained in the oriented intersection graph.
   *
   *   Hence a scc containg more than two nodes is represented as an
   *   array of SccComponent entries,
   *   \f$[ y_{j_0},x_{\phi(j_0)},\dots,y_{j_k},x_{\phi(j_k)}]\f$.
   *
   *   Parameters
   *   scclist ~ resulting sccs
   */

  template<class View, class Tuple>
  forceinline void
  computesccs(Space* home,
	      ViewArray<Tuple>& xz,
	      ViewArray<View>& y,
	      int phi[],
	      SccComponent sinfo[],
	      int scclist[]) {

    // number of sccs is bounded by xs (number of x-nodes)
    int xs = xz.size();
    Support::StaticStack<int> cs(xs);

    //select an y node from the graph
    for (int j = 0; j < xs; j++) {
      int yjmin      = y[j].min();    // the processed min
      while (!cs.empty() && xz[phi[sinfo[cs.top()].rightmost]][0].max() < yjmin) {
	// the topmost scc cannot "reach" y_j or a node to the right of it
	cs.pop();
      }

      // a component has the form C(y-Node, matching x-Node)
      // C is a minimal scc in the oriented intersection graph
      // we only store y_j-Node, since \phi(j) gives the matching X-node
      int i     = phi[j];
      int ximin = xz[i][0].min();
      while (!cs.empty() && ximin <= y[sinfo[cs.top()].rightmost].max()) {
	// y_j can "reach" cs.top() ,
	// i.e. component c can reach component cs.top()
	// merge c and cs.top() into new component
	int top = cs.top();
	// connecting
	sinfo[sinfo[j].leftmost].left        = top;
	sinfo[top].right                     = sinfo[j].leftmost;
	// moving leftmost
	sinfo[j].leftmost                    = sinfo[top].leftmost;
	// moving rightmost
	sinfo[sinfo[top].leftmost].rightmost = j;
	cs.pop();
      }
      cs.push(j);
    }
    cs.reset();


    // now we mark all components with the respective scc-number
    // labeling is bound by O(k) which is bound by O(n)

    for (int i = 0; i < xs; i++) {
      if (sinfo[i].left == i) { // only label variables in sccs
	int scc = sinfo[i].rightmost;
	int z   = i;
	//bound by the size of the largest scc = k
	while (sinfo[z].right != z) {
	  sinfo[z].rightmost   = scc;
	  scclist[phi[z]]      = scc;
	  z                    = sinfo[z].right;
	}
	sinfo[z].rightmost     = scc;
	scclist[phi[z]]        = scc;
      }
    }
  }

  /**
   *   \brief Narrowing the domains of the x variables
   *
   *   Due to the correspondance between perfect matchings in the "reduced"
   *   intersection graph of \a x and \a y views and feasible
   *   assignments for the Sortedness constraint the new domain bounds for
   *   views in \a x are computed as
   *      - lower bounds:
   *        \f$ S_i \geq E_l \f$
   *        where \f$y_l\f$ is the leftmost neighbour of \f$x_i\f$
   *      - upper bounds:
   *        \f$ S_i \leq  E_h \f$
   *        where \f$y_h\f$ is the rightmost neighbour of \f$x_i\f$
   */

  template<class View, class Tuple, bool Perm>
  forceinline bool
  narrow_domx(Space* home,
	      ViewArray<Tuple>& xz,
	      ViewArray<View>& y,
	      int tau[],
	      int phi[],
	      int scclist[],
	      SccComponent sinfo[],
	      bool& nofix) {

    int xs        = xz.size();

    // For every x node
    for (int i = 0; i < xs; i++) {

      int xmin = xz[i][0].min();
      /*
       * take the scc-list for the current x node
       * start from the leftmost reachable y node of the scc
       * and check which Y node in the scc is
       * really the rightmost node intersecting x, i.e.
       * search for the greatest lower bound of x
       */
      int start = sinfo[scclist[i]].leftmost;
      while (y[start].max() < xmin) {
	start = sinfo[start].right;
      }
	
      if (Perm) {
	// start is the leftmost-position for x_i
	// that denotes the lower bound on p_i

	ModEvent me_plb = xz[i][1].gq(home, start);
	if (me_failed(me_plb)) {
	  return false;
	}
	nofix |= (me_modified(me_plb) && start != xz[i][1].min());
      }

      ModEvent me_lb = xz[i][0].gq(home, y[start].min());
      if (me_failed(me_lb)) {
	return false;
      }
      nofix |= (me_modified(me_lb) &&
		y[start].min() != xz[i][0].min());
	
      int ptau = tau[xs - 1 - i];
      int xmax = xz[ptau][0].max();
      /*
       * take the scc-list for the current x node
       * start from the rightmost reachable node and check which
       * y node in the scc is
       * really the rightmost node intersecting x, i.e.
       * search for the smallest upper bound of x
       */
      start = sinfo[scclist[ptau]].rightmost;
      while (y[start].min() > xmax) {
	start = sinfo[start].left;
      }

      if (Perm) {
	//start is the rightmost-position for x_i
	//that denotes the upper bound on p_i
	ModEvent me_pub = xz[ptau][1].lq(home, start);
	if (me_failed(me_pub)) {
	  return false;
	}
	nofix |= (me_modified(me_pub) && start != xz[ptau][1].max());
      }

      ModEvent me_ub = xz[ptau][0].lq(home, y[start].max());
      if (me_failed(me_ub)) {
	return false;
      }
      nofix |= (me_modified(me_ub) &&
		y[start].max() != xz[ptau][0].max());
    }
    return true;
  }

  /**
   * \brief Narrowing the domains of the y views
   *
   * analogously to the x views we take
   *    - for the upper bounds the matching \f$\phi\f$ computed in  glover
   *      and compute the new upper bound by \f$T_j=min(E_j, D_{\phi(j)})\f$
   *    - for the lower bounds the matching \f$\phi'\f$ computed in  revglover
   *      and update the new lower bound by \f$T_j=max(E_j, D_{\phi'(j)})\f$
   */

  template<class View, class Tuple, bool Perm>
  forceinline bool
  narrow_domy(Space* home,
	      ViewArray<Tuple>& xz,
	      ViewArray<View>& y,
	      int phi[],
	      int phiprime[],
	      bool& nofix) {

    int xs = xz.size();
    for (int i = xs; i--; ) {
      ModEvent me_lb = y[i].gq(home, xz[phiprime[i]][0].min());
      if (me_failed(me_lb)) {
	return false;
      }
      nofix |= (me_modified(me_lb) &&
		xz[phiprime[i]][0].min() != y[i].min());

      ModEvent me_ub = y[i].lq(home, xz[phi[i]][0].max());
      if (me_failed(me_ub)) {
	return false;
      }
      nofix |= (me_modified(me_ub) &&
		xz[phi[i]][0].max() != y[i].max());
    }
    return true;
  }

}}}

// STATISTICS: int-prop