/*----------------------------------------------------------------*
 *
 * File : unitClosures.c
 * Author : NTM
 * Created : 12/01/05
 *
 *
 * Copyright (C) Nicolas Thierry-Mieg, 2006.
 *
 *
 * This file is part of InterPool, written by 
 * Nicolas Thierry-Mieg (CNRS, France) Nicolas.Thierry-Mieg@imag.fr
 *
 * InterPool 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 of the License, or
 * (at your option) any later version.
 *
 * InterPool 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.
 *
 * You should have received a copy of the GNU General Public License
 * along with InterPool; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
 *
 *-----------------------------------------------------------------*/

#include <stdio.h> // for printing
#include <stdlib.h> // malloc, exit...
#include <string.h> /* memcpy */

#include "distance.h" /* DIST_XXX*/
#include "types.h" /* MOT */
#include "pools.h" /* pool, getpoolvalue */
#include "closure.h" /* closure */
#include "conflicts.h" /* buildConflictingVars */
#include "poolInVector.h" /* nbMotsPerPoolsVec, setPoolInVector */

#include "unitClosures.h"

/* SANITY: define for extra sanity checks */
#undef SANITY
//#define SANITY

/* DEBUG: define large value for extra output to stderr */
#undef DEBUG
//#define DEBUG 3 


/* BUILD_EMPTY_UNIT_INIT: if defined, initialize memory to 0 in 
   buildEmptyUnitClosures. This wastes a LOT of time (not visible in
   gprof but obvious with time), but could be useful for debugging.
   Definitely undef this in production builds! */
//#define BUILD_EMPTY_UNIT_INIT
#undef BUILD_EMPTY_UNIT_INIT


/* DO_RESIZE_UNITS: if defined, resizeUnitClosures really resizes
   the unitClosures. Otherwise, it only frees them if nbOfClosures==0,
   but leaves the useless memory as is otherwise.
   If you have RAM problems, defining this should help a little.
   Otherwise, undef it so that we don't waste time realloc'ing. */
#undef DO_RESIZE_UNITS

#ifdef DO_RESIZE_UNITS
/* RESIZE_CHECK_NOMOVE: 
   This is meaningless if DO_RESIZE_UNITS is off.
   If defined, check whether memory blocks are moved in resizeUnitClosures. 
   This shouldn't happen (we are just realloc'ing to smaller sizes), 
   and it doesn't on my current system, but it's not illegal.
   If you ever see it happening (bad implementation of realloc?),
   you should definitely undef _DO_RESIZE_UNITS, except if memory is
   a big problem (in which case you should undef RESIZE_CHECK_NOMOVE,
   but be aware that this could cost a LOT of time).
   There is some overhead, so make sure you undef this for production 
   builds once satisfied that your realloc isn't stupid. */
#undef RESIZE_CHECK_NOMOVE
// #define RESIZE_CHECK_NOMOVE
#endif /* DO_RESIZE_UNITS */

/************************************************************************
 ******************* LOCAL FUNCTIONS ************************************
 ************************************************************************/

/////////////// DECLARATIONS ///////////////


/***********************************************************************
 *************** unitC type and functions ******************************
 ***********************************************************************/

/*!
  To build a unitClosures, we first build a set of individual
  unit closures... But the closure datatype doesn't have or
  need faintPools and weakPools which we do need temporarily 
  to build a unitClosures.
  Therefore instead of making the closure type from closure.h
  more complex, we prefer to define the private unitC type, very 
  similar but with the additional required fields.
*/
typedef struct
{
  MOT* negPools ; /*!< conflicting neg pools containing var */
  MOT* faintPools ; /*!< conflicting faint pools containing var */
  MOT* weakPools ; /*!< conflicting weak pools containing var */
  MOT* posPools ; /*!< conflicting pos pools containing var */
  int negPoolsCost ; /*!< total cost of the neg and faint pools in this unitC */
  int score ; /*< distance of all pos and weak pools - negPoolsCost;
		this is (a lower bound of) the score of this closure */
} unitC ;


/*!
  Therefore I must also redefine the set of closures type
*/
typedef struct
{
  int nbOfClosures; /*!< number of closures in the set */
  unitC** allClosuresP; /*!< pointer to the table of unitC pointers. 
			      allClosuresP[i] is a pointer to unitC i, for i < nbOfClosures  */
} setOfUnitCs ;


/*!
  \brief Build and return a pointer to an empty unitC (ie 
  xxxPools are empty (NOT NULL!), and negPoolsCost==score==0).
*/
static unitC* buildEmptyUnitC(int nbPools) ;


/*!
  \brief Free all memory relevant to myUnitC.
*/
static void freeUnitC(unitC* myUnitC) ;


/*!
  \brief Return an empty setOfUnitCs (mem allocated here).<br>
  content is (O, NULL).
*/
static setOfUnitCs* buildEmptySetOfUnitCs(void) ;


/*!
  \brief Free all memory in myUnitCs
*/
static void freeSetOfUnitCs(setOfUnitCs* myUnitCs) ;


/*!
  \brief Given an enriched observation sig, build the unit closure of 
  variable var and return it (mem is allocated here).
*/
static unitC* buildUnitClosure(MOT* tabpool, signature* sig, int var, int n) ;


/*
  Given an enriched observation sig, build each individual unit closure
  (ie the closure of each conflicting variable).
  Return a setOfUnitCs* (nbOfClosures is the number of conflicting vars).
*/
static setOfUnitCs* buildIndivUnitClosures(MOT* tabpool, signature* sig, int n) ;


/*!
  \brief Compare 2 unitC's for use by qsort, to sort the unit closures
  soas to optimize the search speed.
  Currently, sort by decreasing score and then by increasing negPoolsCost.
*/
static int compareUnitCsForSearch(const void* unitC1PP, const void* unitC2PP) ;


/*!
  \brief Given a set of individual unitC's, sort them in order 
  to optimize the search algorithm.
  Currently, we sort by score, and by negPoolsCost when scores 
  are equal, see compareUnitCsForSearch.
*/
static void sortUnitCsForSearch(setOfUnitCs* myClosuresP) ;


/***********************************************************************
 ************ unitClosures functions ***********************************
 ***********************************************************************/


/*!
  \brief Allocate memory for a unitClosures with nbOfClosures closures,
  and initialize nbOfClosures and nbValid (to nbOfClosures). 
  ALL OTHER MEMORY IS ALLOCATED BUT REMAINS UNINITIALIZED, 
  except if you define BUILD_EMPTY_UNIT_INIT above.
  You shouldn't do that, except for debugging.
  Precondition: we must have nbOfClosures > 0  (checked).
*/
static unitClosures* buildEmptyUnitClosures(int nbOfClosures, int nbPools) ;


/*!
  \brief given myUnitsP and a new nbOfClosures (which must 
  be <= myUnitsP->nbOfClosures), "resize" all fields in myUnitsP so
  that it now only contains data for the first nbOfClosures units.
  Actually, resizing only happens if DO_RESIZE_UNITS is defined; otherwise,
  we just leave the useless mem alone (except if nbOfClosures==0). See 
  comments for DO_RESIZE_UNITS and RESIZE_CHECK_NOMOVE above.
  In any case, if nbOfClosures==0 we free everything.
*/
static void resizeUnitClosures(unitClosures* myUnitsP, int nbOfClosures, int nbPools) ;


/*!
  \brief return a pointer to the block of MOTs storing the negPools
  of unit closure unitNum.
  This is used for inclusion-invalidation (in substractNegPools.body.c).
  NOTE: this function doesn't copy the contents! If you modify
  the content of the returned pointer, you will also modify myUnitsP...
  BEWARE!
*/
static MOT* getUnitNegPools(unitClosures* myUnitsP, int unitNum, int nbPools) ;

/*!
  \brief return a pointer to the block of MOTs storing the faintPools
  of unit closure unitNum.
  This is used for inclusion-invalidation (in substractNegPools.body.c).
  NOTE: this function doesn't copy the contents! If you modify
  the content of the returned pointer, you will also modify myUnitsP...
  BEWARE!
*/
static MOT* getUnitFaintPools(unitClosures* myUnitsP, int unitNum, int nbPools) ;


//////////////// BODIES ///////////////////


/***********************************************************************
 *************** unitC type and functions ******************************
 ***********************************************************************/

/*!
  \brief Build and return a pointer to an empty unitC (ie 
  xxxPools are empty (NOT NULL!), and negPoolsCost==score==0).
*/
static unitC* buildEmptyUnitC(int nbPools)
{
  /* allocate all memory */

  unitC* newClosureP = (unitC*)malloc(sizeof(unitC)) ;
  if (newClosureP==NULL)
    {
      fprintf(stderr, "in buildEmptyUnitC, no more memory for newClosureP\n") ;
      exit(1) ;
    }

  int bytesInPools = nbMotsPerPoolsVec(nbPools) * sizeof(MOT) ;

  newClosureP->negPools = (MOT*)malloc(bytesInPools) ;
  if (newClosureP->negPools == NULL)
    {
      fprintf(stderr, "In buildEmptyUnitC, no more memory for negPools\n") ;
      exit(1) ;
    }

  newClosureP->faintPools = (MOT*)malloc(bytesInPools) ;
  if (newClosureP->faintPools == NULL)
    {
      fprintf(stderr, "In buildEmptyUnitC, no more memory for faintPools\n") ;
      exit(1) ;
    }
 
  newClosureP->weakPools = (MOT*)malloc(bytesInPools) ;
  if (newClosureP->weakPools == NULL)
    {
      fprintf(stderr, "In buildEmptyUnitC, no more memory for weakPools\n") ;
      exit(1) ;
    }

  newClosureP->posPools = (MOT*)malloc(bytesInPools) ;
  if (newClosureP->posPools == NULL)
    {
      fprintf(stderr, "In buildEmptyUnitC, no more memory for posPools\n") ;
      exit(1) ;
    }

  /* initialize all memory to 0 */
  memset(newClosureP->negPools, 0, bytesInPools) ;
  memset(newClosureP->faintPools, 0, bytesInPools) ;
  memset(newClosureP->weakPools, 0, bytesInPools) ;
  memset(newClosureP->posPools, 0, bytesInPools) ;
  newClosureP->negPoolsCost = 0 ;
  newClosureP->score = 0 ; 

  return(newClosureP) ;
}

/*!
  \brief Free all memory relevant to myUnitC.
*/
static void freeUnitC(unitC* myUnitC)
{
  /* free Pools */
  free(myUnitC->negPools) ;
  free(myUnitC->faintPools) ;
  free(myUnitC->weakPools) ;
  free(myUnitC->posPools) ;
  /* and free the closure itself */
  free(myUnitC) ;
}


/*!
  \brief Return an empty setOfUnitCs (mem allocated here).<br>
  content is (O, NULL).
*/
static setOfUnitCs* buildEmptySetOfUnitCs(void)
{
  setOfUnitCs* mySetP = (setOfUnitCs*)malloc(sizeof(setOfUnitCs)) ;
  if (mySetP==NULL)
    {
      fprintf(stderr, "in buildEmptySetOfUnitCs, no more memory\n") ;
      exit(1) ;
    }
  mySetP->nbOfClosures = 0 ;
  mySetP->allClosuresP = NULL ;
  return(mySetP) ;
}


/*!
  \brief Free all memory in myUnitCs
*/
static void freeSetOfUnitCs(setOfUnitCs* mySetP)
{
  int i ;
  for (i=0; i < mySetP->nbOfClosures; i++)
    freeUnitC(mySetP->allClosuresP[i]) ;

  free(mySetP->allClosuresP) ;
  free(mySetP) ;
}


/*!
  \brief Given an enriched observation sig, build the unit closure of 
  variable var and return it (mem is allocated here).
*/
static unitC* buildUnitClosure(MOT* tabpool, signature* sig, int var, int n)
{
  int nbPools = sig->nbPools ;
  unitC* myClosureP = buildEmptyUnitC(nbPools) ;

  int absPoolNum ;
  for (absPoolNum=0; absPoolNum < nbPools; absPoolNum++)
    {
      MOT* thispool = pool(tabpool, absPoolNum, n) ;
      if (getpoolvalue(thispool, var) == 11)
	{
	  int sigValue = getsigvalue(sig, absPoolNum) ;
	  /* thispool contains var: see if it is conflicting negative */
	  if  (sigValue == SIG_NEGCONF)
	    {
	      /* it is: include it in negPools, and update negPoolsCost and score */
	      setPoolInVector(myClosureP->negPools, absPoolNum, 1) ;
	      myClosureP->negPoolsCost += DIST_NEG ;
	      myClosureP->score -= DIST_NEG ;
	    }

	  /* else, see if thispool is conflicting faint */
	  else if (sigValue == SIG_FAINTCONF)
	    {
	      setPoolInVector(myClosureP->faintPools, absPoolNum, 1) ;
	      myClosureP->negPoolsCost += DIST_FAINT ;
	      myClosureP->score -= DIST_FAINT ;
	    }
	  
	  /* else, see if thispool is conflicting weak */
	  else if (sigValue == SIG_WEAKCONF)
	    {
	      setPoolInVector(myClosureP->weakPools, absPoolNum, 1) ;
	      myClosureP->score += DIST_WEAK ;
	    }
	  
	  /* else, see if thispool is conflicting positive */
	  else if (sigValue == SIG_POSCONF)
	    {
	      setPoolInVector(myClosureP->posPools, absPoolNum, 1) ;
	      myClosureP->score += DIST_POS ;
	    }
	  
	  /* possible STD optimization: var appears in a single pool in 
	     each layer, therefore we could skip to the next layer.
	     This would just need passing q as a param, or we could
	     calculate it on the fly (largest prime factor of nbPools).
	     Not doing this now, thi function shouldn't be expensive. */
	}
    }
  return(myClosureP) ;
}


/*
  Given an enriched observation sig, build each individual unit closure
  (ie the closure of each conflicting variable).
  Return a setOfUnitCs* (nbOfClosures is the number of conflicting vars).
*/
static setOfUnitCs* buildIndivUnitClosures(MOT* tabpool, signature* sig, int n)
{
  /* to be returned: initially empty set */
  setOfUnitCs* mySetP = buildEmptySetOfUnitCs() ;

  /* allocate enough memory to store all closure pointers no matter how
     many vars are conflicting (ie, n). We will realloc at the end to the
     (smaller) real size */
  mySetP->allClosuresP = (unitC**)malloc(n*sizeof(unitC*)) ;
  if (mySetP->allClosuresP == NULL)
    {
      fprintf(stderr,"in buildIndivUnitClosures, not enough memory for allClosuresP\n");
      exit(1);
    }

  /* identify all conflicting vars. conflictingVars is of size n, 
     element [i] is true if i is a conflicting var. */
  bool* conflictingVars  = buildConflictingVars(tabpool, sig, n) ;


  /* build the unit closure for each variable var */
  int var ;
  for (var=0; var<n; var++)
    {
      /* is var conflicting? if not, skip it */
      if (conflictingVars[var])
	{
	  /* ok, build the unit closure for var */
	  unitC* thisUnitClosureP = buildUnitClosure(tabpool, sig, var, n) ;
	  
	  /* add it to mySetP and next */
	  mySetP->allClosuresP[mySetP->nbOfClosures] = thisUnitClosureP ;
	  mySetP->nbOfClosures++ ;
	}
    }

  /* now resize allClosuresP to just what's needed */
  if (mySetP->nbOfClosures == 0)
    {
      /* deal specifically with this case because return value of realloc with
	 size==0 is not necessarily NULL, but we want to have allClosuresP==NULL */
      free(mySetP->allClosuresP) ;
      mySetP->allClosuresP = NULL ;
    }
  else
    {
      mySetP->allClosuresP = (unitC**)realloc(mySetP->allClosuresP, mySetP->nbOfClosures*sizeof(unitC*)) ;
      if (mySetP->allClosuresP == NULL)
	{
	  fprintf(stderr, "at the end of buildIndivUnitClosures, no more mem for realloc!?\n") ;
	  exit(1) ;
	}
    }

  /* ok, conflicting vars no longer needed, free it */
  free(conflictingVars) ;
  
  return(mySetP);
}


/*!
  Compare 2 unitC's for use by qsort, to sort the unit closures
  soas to optimize the search speed.
  This is for qsort: args are really unitC**'s.
  sort order is defined as follows:
  1. if closure1P->score > closure2P->score, then closure1P < closure2P;
  2. else if the scores are equal:
     2a. if closure1P->negPoolsCost < closure2P->negPoolsCost, closure1P < closure2P;
     2b. else if the negPoolsCost's are equal, then closure1P == closure2P;
  3. otherwise, closure2P < closure1P.

  NOTE 05/04/05: the sort order was changed today (to what is described above), 
  I used to sort by cost first and then by score, but it seems sorting by
  score first is in fact faster. see README.sortOrderForUnitClosures for details.
*/
static int compareUnitCsForSearch(const void* unitC1PP, const void* unitC2PP)
{
  int retval ; /* return value */
  /* convert args to unitC** and get content (which are unitC*) */
  unitC* unitC1P = *(unitC**)unitC1PP ;
  unitC* unitC2P = *(unitC**)unitC2PP ;

  if (unitC1P->score > unitC2P->score)
    retval = -1 ;
  else if (unitC1P->score == unitC2P->score)
    {
      if (unitC1P->negPoolsCost < unitC2P->negPoolsCost)
	retval = -1 ;
      else if (unitC1P->negPoolsCost == unitC2P->negPoolsCost)
	retval = 0;
      else
	retval = 1 ;
    }
  else
    retval = 1 ;

  return(retval) ;
}


/*!
  \brief Given a set of individual unitC's, sort them in order 
  to optimize the search algorithm.
  Currently, we sort by score, and by negPoolsCost when scores 
  are equal, see compareUnitCsForSearch.
*/
static void sortUnitCsForSearch(setOfUnitCs* myClosuresP)
{
  /* simply call qsort, the sort order is defined in compareUnitCsForSearch */
  qsort(myClosuresP->allClosuresP, 
	myClosuresP->nbOfClosures, 
	sizeof(unitC*), compareUnitCsForSearch) ;
}



/***********************************************************************
 ************ unitClosures functions ***********************************
 ***********************************************************************/

/*!
  Allocate memory for a unitClosures with nbOfClosures closures,
  and initialize nbOfClosures and nbValid (to nbOfClosures). 
  ALL OTHER MEMORY IS ALLOCATED BUT REMAINS UNINITIALIZED, 
  except if you define BUILD_EMPTY_UNIT_INIT above.
  You shouldn't do that, except for debugging.
  Precondition: we must have nbOfClosures > 0 (checked).
  This is because malloc(0) is implementation-defined (can return NULL).
*/
static unitClosures* buildEmptyUnitClosures(int nbOfClosures, int nbPools)
{
  /* check precondition */
  if (nbOfClosures <= 0)
    {
      fprintf(stderr, 
	      "buildEmptyUnitClosures called with nbOfClosures (%d) <= 0, you shouldn't need that\n",
	      nbOfClosures) ;
      exit(1) ;
    }

  unitClosures* myClosuresP = (unitClosures*)malloc(sizeof(unitClosures)) ;
  if (myClosuresP==NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for myClosures\n") ;
      exit(1) ;
    }

  myClosuresP->nbOfClosures = nbOfClosures ;
  myClosuresP->nbValid = nbOfClosures ;
  
  /* number of bytes in xxxPools */
  int poolsSize = nbOfClosures * nbMotsPerPoolsVec(nbPools) * sizeof(MOT) ;
  /* number of bytes in the int vectors */
  int intVectorSize = nbOfClosures * sizeof(int) ;

  /* Memory allocation for negPools */
  myClosuresP->negPools = (MOT*)malloc(poolsSize) ;
  if (myClosuresP->negPools == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for negPools\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->negPools, 0, poolsSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */
  

  /* Memory allocation for faintPools */
  myClosuresP->faintPools = (MOT*)malloc(poolsSize) ;
  if (myClosuresP->faintPools == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for faintPools\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->faintPools, 0, poolsSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */
  
  /* Memory allocation for negPoolCosts */
  myClosuresP->negPoolCosts = (int*)malloc(intVectorSize) ;
  if (myClosuresP->negPoolCosts == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for negPoolCosts\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->negPoolCosts, 0, intVectorSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */

  /* Memory allocation for weakPools */
  myClosuresP->weakPools = (MOT*)malloc(poolsSize) ;
  if (myClosuresP->weakPools == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for weakPools\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->weakPools, 0, poolsSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */

  /* Memory allocation for posPools */
  myClosuresP->posPools = (MOT*)malloc(poolsSize) ;
  if (myClosuresP->posPools == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for posPools\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->posPools, 0, poolsSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */

  /* Memory allocation for scores */
  myClosuresP->scores = (int*)malloc(intVectorSize) ;
  if (myClosuresP->scores == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for scores\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->scores, 0, intVectorSize) ;
#endif /* BUILD_EMPTY_UNIT_INIT */

  /* Memory allocation for valid */
  myClosuresP->valid = (bool*)malloc(nbOfClosures * sizeof(bool)) ;
  if (myClosuresP->valid == NULL)
    {
      fprintf(stderr, "in buildEmptyUnitClosures, no more memory for valid\n") ;
      exit(1) ;
    }
#ifdef BUILD_EMPTY_UNIT_INIT
  memset(myClosuresP->valid, 0, nbOfClosures * sizeof(bool)) ;
#endif /* BUILD_EMPTY_UNIT_INIT */

  return(myClosuresP) ;
}


/*
  \brief given myUnitsP and a new nbOfClosures (which must 
  be <= myUnitsP->nbOfClosures), "resize" all fields in myUnitsP so
  that it now only contains data for the first nbOfClosures units.
  Actually, resizing only happens if DO_RESIZE_UNITS is defined; otherwise,
  we just leave the useless mem alone (except if nbOfClosures==0). See 
  comments for DO_RESIZE_UNITS and RESIZE_CHECK_NOMOVE above.
  In any case, if nbOfClosures==0 we free everything.

  Also, in DO_RESIZE_UNITS mode we use realloc; realloc specifies that the 
  content doesn't change, but the pointers can change (in this case, the 
  content is copied... slowly). In our case, this shouldn't happen because 
  we always realloc to a smaller size... But if the realloc implementation 
  is bad it could happen! To test for this, define RESIZE_CHECK_NOMOVE.
*/
static void resizeUnitClosures(unitClosures* myUnitsP, int nbOfClosures, int nbPools)
{
#ifdef SANITY
  if (myUnitsP->nbOfClosures < nbOfClosures)
    {
      fprintf(stderr, 
	      "in resizeUnitClosures, myUnitsP->nbOfClosures (==%d) < nbOfClosures (==%d)!\n",
	      myUnitsP->nbOfClosures, nbOfClosures) ;
      exit(1) ;
    }
#endif /* SANITY */

#ifndef DO_RESIZE_UNITS
  /* Don't resize at all.
     We still have to free memory if nbOfClosures==0, because
     freeUnitClosures expects that.
     But in general, this should improve performance (realloc'ing
     costs something), at the cost of some RAM (the useless memory
     is not freed straight away).
     Since RAM consumption is currently low, that shouldn't matter
     too much, but in case of problems define DO_RESIZE_UNITS.
  */
  if ((nbOfClosures==0) && (myUnitsP->nbOfClosures != 0))
    {
      /* free the memory, and set pointers to NULL for safety */
      free(myUnitsP->negPools) ;
      myUnitsP->negPools = NULL ;
      free(myUnitsP->faintPools) ;
      myUnitsP->faintPools = NULL ;
      free(myUnitsP->negPoolCosts) ;
      myUnitsP->negPoolCosts = NULL ;
      free(myUnitsP->weakPools) ;
      myUnitsP->weakPools = NULL ;
      free(myUnitsP->posPools) ;
      myUnitsP->posPools = NULL ;
      free(myUnitsP->scores) ;
      myUnitsP->scores = NULL ;
      free(myUnitsP->valid) ;
      myUnitsP->valid = NULL ;
    }

  /* in any case, set myUnitsP->nbOfClosures */
  myUnitsP->nbOfClosures = nbOfClosures ;
  /* nbValid shouldn't change */

  /* silence compiler warnings on unused nbPools */
  if (FALSE)
    {
      int bogus = nbPools ;
      bogus++ ;
    }
  
#else /* DO_RESIZE_UNITS is on */

#ifdef RESIZE_CHECK_NOMOVE
  /* save old pointers, to see if memory gets moved */
  MOT* oldNegPools = myUnitsP->negPools ;
  MOT* oldFaintPools = myUnitsP->faintPools ;
  int* oldNegCosts = myUnitsP->negPoolCosts ;
  MOT* oldWeakPools = myUnitsP->weakPools ;
  MOT* oldPosPools = myUnitsP->posPools ;
  int* oldScores = myUnitsP->scores ;
  bool* oldValid = myUnitsP->valid ;
#endif /* RESIZE_CHECK_NOMOVE */

  /* new number of bytes in xxxPools */
  int newPoolsSize = nbOfClosures * nbMotsPerPoolsVec(nbPools) * sizeof(MOT) ;
  /* new number of bytes in the int vectors */
  int newIntVectorSize = nbOfClosures * sizeof(int) ;

  /* do reallocations */
  myUnitsP->nbOfClosures = nbOfClosures ;
  myUnitsP->negPools = realloc(myUnitsP->negPools, newPoolsSize) ;
  if (myUnitsP->negPools == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc negPools\n") ;
      exit(1) ;
    }
  myUnitsP->faintPools = realloc(myUnitsP->faintPools, newPoolsSize) ;
  if (myUnitsP->faintPools == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc faintPools\n") ;
      exit(1) ;
    }
  myUnitsP->negPoolCosts = realloc(myUnitsP->negPoolCosts, newIntVectorSize) ;
  if (myUnitsP->negPoolCosts == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc negPoolCosts\n") ;
      exit(1) ;
    }
  myUnitsP->weakPools = realloc(myUnitsP->weakPools, newPoolsSize) ;
  if (myUnitsP->weakPools == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc weakPools\n") ;
      exit(1) ;
    }
  myUnitsP->posPools = realloc(myUnitsP->posPools, newPoolsSize) ;
  if (myUnitsP->posPools == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc posPools\n") ;
      exit(1) ;
    }
  myUnitsP->scores = realloc(myUnitsP->scores, newIntVectorSize) ;
  if (myUnitsP->scores == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc scores\n") ;
      exit(1) ;
    }
  myUnitsP->valid = realloc(myUnitsP->valid, nbOfClosures * sizeof(bool)) ;
  if (myUnitsP->negPools == NULL)
    {
      fprintf(stderr, "in resizeUnitClosures, no more mem for realloc valid\n") ;
      exit(1) ;
    }

  /* if nbOfClosures == 0: set pointers to NULL (realloc doesn't have to do it) */
  if (nbOfClosures==0)
    {
      myUnitsP->negPools = NULL ;
      myUnitsP->faintPools = NULL ;
      myUnitsP->negPoolCosts = NULL ;
      myUnitsP->weakPools = NULL ;
      myUnitsP->posPools = NULL ;
      myUnitsP->scores = NULL ;
      myUnitsP->valid = NULL ;
    }
#ifdef RESIZE_CHECK_NOMOVE
  /* otherwise, test whether pointers changed (== content was moved) */
  else
    {
      bool contentMoved = FALSE ;
      if (oldNegPools != myUnitsP->negPools)
	{
	  fprintf(stderr, "in resizeUnitClosures, negPools got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldFaintPools != myUnitsP->faintPools)
	{
	  fprintf(stderr, "in resizeUnitClosures, faintPools got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldNegCosts != myUnitsP->negPoolCosts)
	{
	  fprintf(stderr, "in resizeUnitClosures, negPoolCosts got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldWeakPools != myUnitsP->weakPools)
	{
	  fprintf(stderr, "in resizeUnitClosures, weakPools got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldPosPools != myUnitsP->posPools)
	{
	  fprintf(stderr, "in resizeUnitClosures, posPools got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldScores != myUnitsP->scores)
	{
	  fprintf(stderr, "in resizeUnitClosures, scores got moved!\n") ;
	  contentMoved = TRUE ;
	}
      if (oldValid != myUnitsP->valid)
	{
	  fprintf(stderr, "in resizeUnitClosures, valid got moved!\n") ;
	  contentMoved = TRUE ;
	}

      if (contentMoved)
	exit(1) ;
    }
#endif /* RESIZE_CHECK_NOMOVE */

#endif /* DO_RESIZE_UNITS */
}


/*
  \brief return a pointer to the block of MOTs storing the negPools
  of unit closure unitNum.
  This is used for inclusion-invalidation (in substractNegPools.body.c).
  NOTE: this function doesn't copy the contents! If you modify
  the content of the returned pointer, you will also modify myUnitsP...
  BEWARE!
*/
static MOT* getUnitNegPools(unitClosures* myUnitsP, int unitNum, int nbPools)
{
  return(myUnitsP->negPools + (unitNum * nbMotsPerPoolsVec(nbPools))) ;
}

/*
  Return a pointer to the block of MOTs storing the faintPools
  of unit closure unitNum.
  This is used for inclusion-invalidation (in substractNegPools.body.c).
  NOTE: this function doesn't copy the contents! If you modify
  the content of the returned pointer, you will also modify myUnitsP...
  BEWARE!
*/
static MOT* getUnitFaintPools(unitClosures* myUnitsP, int unitNum, int nbPools)
{
  return(myUnitsP->faintPools + (unitNum * nbMotsPerPoolsVec(nbPools))) ;
}


/************************************************************************
 ******************* EXPORTED FUNCTIONS *********************************
 ************************************************************************/

/*
  free all memory from myClosuresP
*/
void freeUnitClosures(unitClosures* myClosuresP)
{
  if (myClosuresP->nbOfClosures != 0)
    {
      /* free memory; if nbOfClosures==0, nothing to free */
      free(myClosuresP->negPools) ;
      free(myClosuresP->faintPools) ;
      free(myClosuresP->negPoolCosts) ;
      free(myClosuresP->weakPools) ;
      free(myClosuresP->posPools) ;
      free(myClosuresP->scores) ;
      free(myClosuresP->valid) ;
    }
  free(myClosuresP) ;
}

/*
  Return a copy of myClosuresP. mem is allocated here.
*/
unitClosures* copyUnitClosures(unitClosures* myClosuresP, int nbPools)
{
  unitClosures* newUnitsP = (unitClosures*)malloc(sizeof(unitClosures));

  int nbOfClosures = myClosuresP->nbOfClosures;
  /* number of bytes in xxxPools */
  int poolsSize = nbOfClosures * nbMotsPerPoolsVec(nbPools) * sizeof(MOT) ;
  /* number of bytes in the int vectors */
  int intVectorSize = nbOfClosures * sizeof(int) ;

  newUnitsP->nbOfClosures = nbOfClosures ;
  newUnitsP->nbValid = myClosuresP->nbValid ;

  /* Memory allocation for negPools*/
  newUnitsP->negPools = (MOT*)malloc(poolsSize) ;
  if (newUnitsP->negPools == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for negPools\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->negPools, myClosuresP->negPools, poolsSize) ;
  
  /* Memory allocation for faintPools*/
  newUnitsP->faintPools = (MOT*)malloc(poolsSize) ;
  if (newUnitsP->faintPools == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for faintPools\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->faintPools, myClosuresP->faintPools, poolsSize) ;
  
  /* Memory allocation for negPoolCosts*/
  newUnitsP->negPoolCosts = (int*)malloc(intVectorSize) ;
  if (newUnitsP->negPoolCosts == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for negPoolCosts\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->negPoolCosts, myClosuresP->negPoolCosts, intVectorSize) ;

  /* Memory allocation for weakPools */
  newUnitsP->weakPools = (MOT*)malloc(poolsSize) ;
  if (newUnitsP->weakPools == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for weakPools\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->weakPools, myClosuresP->weakPools, poolsSize) ;

  /* Memory allocation for posPools */
  newUnitsP->posPools = (MOT*)malloc(poolsSize) ;
  if (newUnitsP->posPools == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for posPools\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->posPools, myClosuresP->posPools, poolsSize) ;

  /* Memory allocation for scores */
  newUnitsP->scores = (int*)malloc(intVectorSize) ;
  if (newUnitsP->scores == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for scores\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->scores, myClosuresP->scores, intVectorSize) ;


  /* Memory allocation for valid */
  newUnitsP->valid = (bool*)malloc(nbOfClosures * sizeof(bool)) ;
  if (newUnitsP->valid == NULL)
    {
      fprintf(stderr, "in copyUnitClosures, no more memory for valid\n") ;
      exit(1) ;
    }
  memcpy(newUnitsP->valid, myClosuresP->valid, nbOfClosures * sizeof(bool)) ;


  return(newUnitsP);
}


/*
  Return a copy of myUnitsP->valid. mem is allocated here.
*/
bool* copyValidVector(unitClosures* myUnitsP)
{
  int bytesInValid = myUnitsP->nbOfClosures * sizeof(bool) ;
  /* allocate a new vector */
  bool* newValidVector = (bool*)malloc(bytesInValid) ;
  if (newValidVector == NULL)
    {
      fprintf(stderr, "in copyValidVector, no more memory for new valid vector\n") ;
      exit(1) ;
    }
  /* copy content of myUnits->valid */
  memcpy(newValidVector, myUnitsP->valid, bytesInValid) ;
  return(newValidVector) ;
}


/*!
  \brief Perform expensive checks on myUnitsP:
  - Compare nbValid with the number of TRUE's in valid;
  - Count the number of non-zero bits in each unit's negPools
  and faintPools, and multiply by corresponding DIST_xxx: 
  this should be equal to the unit's negPoolCosts.
  Die with error message if a discrepancy is found.
  This is for sanity-checking.
  NOTE: this function relies on the fact that xxxPools are
  padded with 0's
*/
void checkUnits(unitClosures* myUnitsP, int nbPools)
{
  /* first, test nbValid */
  {
    bool* validVector = myUnitsP->valid ;
    int numValid = 0 ; /* number of valid units */
    int i ;
    for (i=myUnitsP->nbOfClosures; i>0; i--)
      {
	if (*validVector)
	  numValid++ ;
	validVector++ ;
      }
    if (numValid != myUnitsP->nbValid)
      {
	fprintf(stderr, 
		"in checkUnits, numValid (%d) <= nbValid (%d)!\n",
		numValid, myUnitsP->nbValid) ;
	exit(1) ;
      }
  }
  
  /* now test negPoolCosts */
  {
    /* number of MOTs to store a pools vector: */
    int motsInPools = nbMotsPerPoolsVec(nbPools) ;

    MOT* negPools = myUnitsP->negPools ;
    MOT* faintPools = myUnitsP->faintPools ;
    int* costs = myUnitsP->negPoolCosts ;

    int unitCounter ;
    for (unitCounter = 0; unitCounter < myUnitsP->nbOfClosures; unitCounter++)
      {
	int calcCost = 0 ; /* calculated cost */
	int motCounter ;
	for (motCounter = motsInPools ; motCounter > 0 ; motCounter--)
	  {
	    MOT currentMot = *negPools ;
	    while (currentMot != 0)
	      {
		/* clear the lowest weight bit at 1 in currentMot */
		currentMot &= (currentMot-1) ;
		calcCost += DIST_NEG ;
	      }
	    /* same for faint */
	    currentMot = *faintPools ;
	    while (currentMot != 0)
	      {
		/* clear the lowest weight bit at 1 in currentMot */
		currentMot &= (currentMot-1) ;
		calcCost += DIST_FAINT ;
	      }
	    negPools++ ;
	    faintPools++ ;
	  }
	
	if (calcCost != *costs)
	  {
	    fprintf(stderr,
		    "in checkUnits: for unit %d, calculated cost is %d, should be %d\n",
		    unitCounter, calcCost, *costs) ;
	    exit(1) ;
	  }
	
	costs++ ;
      }
  }
}


/*!
  \brief Given an enriched observation sig, build all the unit closures
  (ie the closures of each single conflicting variable), sort them 
  in an a priori favorable order (as specified in sortClosuresForSearch), 
  and return them in a (ptr to a) unitClosures.
  pre-condition: there must be at least one conflicting var in sig.

  Algorithm is:
  - 1. build individual unit closures
  - 2. sort them appropriately
  - 3. copy the contents into a unitClosures structure
  - 4. free the individual closures and return
*/
unitClosures* buildUnitClosures(MOT* tabpool, signature* sig, int n)
{
  int nbPools = sig->nbPools ;


  /* 1. build individual closures */
  setOfUnitCs* indivUnitClosuresP = buildIndivUnitClosures(tabpool, sig, n) ;

  /* shouldn't call this function if there aren't any conflicting
     vars, but just in case: */
  if (indivUnitClosuresP->nbOfClosures == 0)
    {
      fprintf(stderr, "buildUnitClosures called but no conflicts in sig\n") ;
      exit(1) ;
    }

  /* 2. sort the closures */
  sortUnitCsForSearch(indivUnitClosuresP) ;


  /* 3. copy content to myClosuresP which will be returned */
  /* allocate mem and set nbOfClosures and nbValid */
  unitClosures* myClosuresP = buildEmptyUnitClosures(indivUnitClosuresP->nbOfClosures, nbPools) ;
  
  /* copy all fields of each unitC, and set valid  */
  {
    /* maintain pointers to the current space to fill in each table */
    MOT* currentNeg = myClosuresP->negPools ;
    MOT* currentFaint = myClosuresP->faintPools ;
    int* currentNegCost = myClosuresP->negPoolCosts ;
    MOT* currentWeak = myClosuresP->weakPools ;
    MOT* currentPos = myClosuresP->posPools ;
    int* currentScore = myClosuresP->scores ;
    bool* currentValid = myClosuresP->valid ;

    int blockSize = nbMotsPerPoolsVec(nbPools) ; /* number of MOTs in a block */
    int blockBytes = blockSize * sizeof(MOT) ; /* number of bytes in a block */

    int var ;
    for (var=0; var < indivUnitClosuresP->nbOfClosures; var++)
      {
	unitC* thisUnitP = indivUnitClosuresP->allClosuresP[var];
	
	memcpy(currentNeg, thisUnitP->negPools, blockBytes) ;
	currentNeg += blockSize ;

	memcpy(currentFaint, thisUnitP->faintPools, blockBytes) ;
	currentFaint += blockSize ;
	
	*currentNegCost = thisUnitP->negPoolsCost ;
	currentNegCost++ ;
		
	memcpy(currentWeak,thisUnitP->weakPools, blockBytes) ;
	currentWeak += blockSize ;
	
	memcpy(currentPos,thisUnitP->posPools, blockBytes) ;
	currentPos += blockSize ;
	
	*currentScore = thisUnitP->score ;
	currentScore++ ;
	
	*currentValid = TRUE ;
	currentValid++ ;
      }
  }

  /* 4. free individual closures */
  freeSetOfUnitCs(indivUnitClosuresP) ;

  return(myClosuresP) ;
}


/*
  copy sourceClosureP's negPools, negPoolsCost and score into 
  destClosureP,squashing previous content.
  Then update these 3 fields from destClosureP by doing the union with
  the negPools and faintPools for unit vector unitNum.
  NOTE: after the call, destClosureP is no longer coherent: it's posPools
  contains garbage, and destClosureP->score is in reality:
  <cost of sourceClosureP->posPools> - destClosureP->negPoolsCost.
  You must either discard *destClosureP, or call finishUnionWithUnit using
  the same params.
  This is somewhat ugly, but it allows to examine destClosureP->negPoolsCost
  for potential pruning before copying and doing the union with the unit's
  posPools and weakPools...
*/
void beginUnionWithUnit(closure* destClosureP, closure* sourceClosureP, 
			unitClosures* myUnitsP, int unitNum, int nbPools)
{
  /* number of MOTs to store a pools vector: */
  int numMotsInPools = nbMotsPerPoolsVec(nbPools) ;
  /* copy sourceClosure's negPoolsCost and score */
  destClosureP->negPoolsCost = sourceClosureP->negPoolsCost ;
  destClosureP->score = sourceClosureP->score ;
  
  /* do the union with the unit's negpools and faintpools */
  MOT* negUnitPools = (myUnitsP->negPools) + (unitNum * numMotsInPools) ;
  MOT* faintUnitPools = (myUnitsP->faintPools) + (unitNum * numMotsInPools) ;
  MOT* sourcePools = sourceClosureP->negPools ;
  MOT* destPools = destClosureP->negPools ;

  int motCounter ;
  for (motCounter = numMotsInPools ; motCounter > 0 ; motCounter--)
    {
      /* do the union with neg pools */
      *destPools = *sourcePools | *negUnitPools ;
      /* update negPoolsCost and score, ie count the number of bits that 
	 changed between src and dst: 
	diffBits contains 1's wherever a pool was added to src */
      MOT diffBits = ~(*sourcePools) & *negUnitPools ;
      while(diffBits != 0)
	{
	  /* clear the lowest-weight 1 bit in diffBits */
	  diffBits &= (diffBits-1) ;
	  destClosureP->negPoolsCost += DIST_NEG ;
	  destClosureP->score -= DIST_NEG ;
	}
      
      /* do the union with faint pools */
      *destPools = *destPools | *faintUnitPools ;
      /* again, update negPoolsCost and score, but this time the
	 cost of each bit change is DIST_FAINT */
      diffBits = ~(*sourcePools) & *faintUnitPools ;
      while(diffBits != 0)
	{
	  /* clear the lowest-weight 1 bit in diffBits */
	  diffBits &= (diffBits-1) ;
	  destClosureP->negPoolsCost += DIST_FAINT ;
	  destClosureP->score -= DIST_FAINT ;
	}
      
      negUnitPools++ ;
      faintUnitPools++ ;
      sourcePools++ ;
      destPools++ ;
    }
}


/*
  Finish filling destClosureP with the union of sourceClosureP 
  with unit closureunitNum: copy sourceClosureP's posPools into 
  destClosureP, do the union with the unit closure's posPools and
  weakPools, and update score accordingly.
  This function should only be called after beginUnionWithUnit, 
  on the same args.
*/
void finishUnionWithUnit(closure* destClosureP, closure* sourceClosureP, 
			 unitClosures* myUnitsP, int unitNum, int nbPools)
{
  /* number of MOTs to store a pools vector: */
  int numMotsInPools = nbMotsPerPoolsVec(nbPools) ;

  /* do the union of source.pospools with the unit's pospools and
     weakpools, place it into dest.pospools, and update score */
  MOT* weakUnitPools = (myUnitsP->weakPools) + (unitNum * numMotsInPools) ;
  MOT* posUnitPools = (myUnitsP->posPools) + (unitNum * numMotsInPools) ;
  MOT* sourcePools = sourceClosureP->posPools ;
  MOT* destPools = destClosureP->posPools ;

  int motCounter ;
  for (motCounter = numMotsInPools ; motCounter > 0 ; motCounter--)
    {
      /* do the union with weak pools, and update score */
      *destPools = *sourcePools | *weakUnitPools ;
      MOT diffBits = ~(*sourcePools) & *weakUnitPools ;
      while(diffBits != 0)
	{
	  /* clear the lowest-weight 1 bit in diffBits */
	  diffBits &= (diffBits-1) ;
	  destClosureP->score += DIST_WEAK ;
	}
      
      /* do the union with pos pools, and update score */
      *destPools = *destPools | *posUnitPools ;
      diffBits = ~(*sourcePools) & *posUnitPools ;
      while(diffBits != 0)
	{
	  /* clear the lowest-weight 1 bit in diffBits */
	  diffBits &= (diffBits-1) ;
	  destClosureP->score += DIST_POS ;
	}
      
      posUnitPools++ ;
      weakUnitPools++ ;
      sourcePools++ ;
      destPools++ ;
    }
}


/* bodies of substractNegPoolsFromUnits and substractNegPoolsFromUnitsInclInval
   are in substractNegPools.body.
   This is done because these 2 functions are almost identical, but I
   want to have both functions available because inclusion-invalidation
   has a substantial overhead and is rarely required. I see more than 10%
   time wasted when the func is compiled with inclusion-invalidation,
   even when that feature is never used (ie when INCL_INVAL_MINUNITS=100000,
   but you can't see that anymore because the code has changed).
   A script takes care of appending the body of substractNegPoolsFromUnitsInclInval
   at the end of substractNegPools.body, it is called from the Makefile.
*/
#include "substractNegPools.bodyFull.c"


/*!
  \brief Substract currentClosureP->posPools from each valid unit's posPools.
  If any unit's posPools becomes empty, invalidate it.
  This function should be called after substractNegPoolsFromUnits, on the
  unitClosures produced by it.
  If myUnitsP->nbOfClosures==0, it's OK, nothing happens.

  NOTE: we don't bother checking that the units are initially valid anymore: we 
  just do all the substractions, and invalidate the units if possible. We do this
  because this function is called on the unitClosures produced by
  substractNegPoolsFromUnits, where all units are valid (any unit which is
  or becomes invalid, is trashed).
*/
void substractPosPoolsFromUnits(closure* currentClosureP, unitClosures* myUnitsP, int nbPools)
{
  int nbOfClosures = myUnitsP->nbOfClosures ;

  MOT* unitPosPools = myUnitsP->posPools ;
  MOT* unitWeakPools = myUnitsP->weakPools ;
  int* unitScores = myUnitsP->scores ;
  bool* unitValid = myUnitsP->valid ;
  int* unitNegPoolCosts = myUnitsP->negPoolCosts ; /* for checking if pos+weak still non-empty */


#ifdef DEBUG
  int nbInval = 0 ; /* number of units invalidated */
#endif /* DEBUG */

  int unitCounter ; /* counter for the unit closure being updated */
  int motsInPools = nbMotsPerPoolsVec(nbPools) ; /* number of MOTs to store a pools vector */

  for (unitCounter=nbOfClosures; unitCounter>0; unitCounter--)
    {
      
      /* substract currentClosureP->posPools from this unit's posPools and weakPools */
      MOT* currentPosPools = currentClosureP->posPools ;
 
      int motCounter ;
      for (motCounter = motsInPools ; motCounter > 0 ; motCounter--)
	{
	  /* calculate diffBits to update score (BEFORE modifying unitPosPools!) */
	  MOT diffBits = (*unitPosPools) & (*currentPosPools) ;
	  /* now substract currentPosPools from unitPos */
	  *unitPosPools = (*unitPosPools) & ~(*currentPosPools) ;
	  /* update score */
	  while(diffBits != 0)
	    {
	      /* clear the lowest-weight 1 bit in diffBits */
	      diffBits &= (diffBits-1) ;
	      (*unitScores) -= DIST_POS ;
	    }
	  
	  /* same thing for unitWeakPools */
	  diffBits = (*unitWeakPools) & (*currentPosPools) ;
	  *unitWeakPools = (*unitWeakPools) & ~(*currentPosPools) ;
	  while(diffBits != 0)
	    {
	      /* clear the lowest-weight 1 bit in diffBits */
	      diffBits &= (diffBits-1) ;
	      (*unitScores) -= DIST_WEAK ;
	    }

	  /* next MOT */
	  currentPosPools++ ;
	  unitPosPools++ ;
	  unitWeakPools++ ;
	}

      /* if this unit's posPools has become empty, invalidate the unit */
      if ((*unitScores) + (*unitNegPoolCosts) == 0)
	{
	  *unitValid = FALSE ;
	  myUnitsP->nbValid-- ;
#ifdef DEBUG
	  nbInval++ ;
#endif /* DEBUG */
	}
      	

      unitScores++ ;
      unitValid++ ;
      unitNegPoolCosts++ ;
    }

#ifdef DEBUG
  fprintf(stderr, "in substractPosPoolsFromUnits: invalidated %d in new\n", nbInval) ;
#endif /* DEBUG */
}

