eg_greedysample.c

Go to the documentation of this file.

#include "eg_greedysample.h"

/* CHANGES TO EG_GREEDYKP:

   1. Add eo_dist to EGgreedyData_t
   2. Add an array: double *cand_weight, of size (cycleG->nedges).
   3. Add an array: EGdGraphEdge_t **cand_edge, of size (cycleG->nedges).

*/

int EGgenerateRandomInternalPairs( EGgreedyData_t *const gdata,
                                   EGlist_t *const new_pairs,
                                   EGlist_t *const old_pairs,
                                   unsigned int max_new_pairs,
                                   double sample_time )
{

  int rval;
  EGlistNode_t *it;
  EGinternalPairs_t *old_p, *new_p;
  EGdijkstraCost_t ext_ubound;
  EGlist_t *ext_path;
  EGdijkstraCost_t ext_path_val;
  EGmemPool_t *mem = gdata->bdG->mem;
  EGtimer_t timer;

  const int ext_parity = 0;
  const EGdijkstraCost_t dij_two = EGdijkstraToCost(2.0);

  for(it=old_pairs->end; it; it=it->prev)
  {

    old_p = (EGinternalPairs_t*)(it->this);
    ext_ubound = EGdijkstraCostSub(dij_two,old_p->ivalue);
    EGtimerReset(&timer);

    while(1)
    {

      EGtimerStart(&timer);

      ext_path = EGnewList(mem);  
 
      rval = EGgreedySample(old_p->s,
                            old_p->t,
                            ext_ubound,
                            ext_parity,
                            gdata,
                            &ext_path_val,
                            ext_path);
      CHECKRVAL(rval);

      if (!ext_path->size)
      {
        //fprintf(stderr, "ignoring.\n");
        EGfreeList(ext_path);
      }
      else
      {
        //fprintf(stderr, "opt = %lf  new = %lf\n", EGdijkstraCostToLf(old_p->evalue), EGdijkstraCostToLf(ext_path_val));
        new_p = EGnewRandomInternalPairs(old_p, ext_path, ext_path_val, mem);
        TEST(!new_p, "bad internal pair generated");
        EGaddIPtoList(new_p, new_pairs, max_new_pairs); 
      }

      EGtimerStop(&timer);
      //fprintf(stderr, "%lf / %lf\n", timer.time, (sample_time / old_pairs->size) );
      if (timer.time > (sample_time / old_pairs->size) )
         break;
    }

  }

  return 0;

}

EGinternalPairs_t* EGnewRandomInternalPairs( EGinternalPairs_t *old_pair,
                                             EGlist_t *ext_path,
                                             EGdijkstraCost_t ext_value,
                                             EGmemPool_t *mem )
{
  EGinternalPairs_t *p;

  p = EGmemPoolSMalloc(mem, EGinternalPairs_t, 1);

  p->s = old_pair->s;
  p->t = old_pair->t;

  p->value  = old_pair->ivalue + ext_value;
  p->ivalue = old_pair->ivalue;
  p->evalue = ext_value;

  p->npath = old_pair->npath;
  p->stpath = EGmemPoolSMalloc(mem, EGdGraphEdge_t*, p->npath);
  memcpy(p->stpath, old_pair->stpath, sizeof(EGdGraphEdge_t*)*(p->npath));

  p->extpath = ext_path;

  return p;
}

/* NOTE: s,t are given from the graph bdG. 
         we are sampling in the graph cycleG. */

int EGgreedySample( EGdGraphNode_t *bdg_s,
                    EGdGraphNode_t *bdg_t,
                    EGdijkstraCost_t st_ubound, 
                    int st_parity,
                    EGgreedyData_t *gdata,
                    EGdijkstraCost_t *path_val,
                    EGlist_t *path )
{

  unsigned int cnt = 0;
  EGdGraphNode_t *s, *t;
  EGdGraphNode_t *curr_s;
  EGdGraphEdge_t *curr_e;
  EGcycleData_t *cd;
  int curr_parity = 0, rem_parity=st_parity;
  int e_parity = 0;
  int curr_s_bdg_id;
  EGdijkstraCost_t e_val, rem_ubound, curr_val;
  EGlist_t *curr_path;
   
  if (!st_parity)
  {
    s = gdata->bdGtoCycleGmap[ bdg_s->id ];
    t = gdata->bdGtoCycleGmap[ bdg_t->id ];
  }
  else
  {
    s = gdata->bdGtoCycleGmap[ bdg_s->id ];
    t = (EGdGraphNode_t*) gdata->bdGtoCycleGmap[ bdg_t->id ]->cn->next->this;
  }

  curr_path = EGnewList(gdata->bdG->mem);
  EGlistClear(path, nullFree);

  *path_val = st_ubound;

  // here we could iterate if we wanted to, by adding a while or for.
  {

    //fprintf(stderr, "\n*** s = %d, t = %d *** \n\n", s->id, t->id);

    /* initialize */
    curr_s = s; 
    curr_val = EGdijkstraToCost(0.0);
    curr_parity = 0;
    rem_ubound = st_ubound;
    rem_parity = st_parity;
    EGlistClear(curr_path,nullFree);
    memset(gdata->bdnodes, 0, sizeof(EGdGraphNode_t*)*gdata->bdG->nodes->size);

    /* iterate through edges */
    for(cnt=0; cnt < gdata->cycleG->nodes->size; cnt++)
    {

      /* mark the current node as visited */
      curr_s_bdg_id = gdata->bdGnodeMap[curr_s->id/2]->id;
      gdata->bdnodes[curr_s_bdg_id] = (EGdGraphNode_t*)(1);

      /* sample an edge */
      curr_e = EGgreedySampleEdge(curr_s,  
                                  t, 
                                  rem_ubound, 
                                  rem_parity, 
                                  gdata);
      if (!curr_e)
      {
        EGlistClear(curr_path, nullFree);
        //fprintf(stderr, "[dead end]");
        break;
      }

      /* update current path information */
      cd = (EGcycleData_t*)curr_e->data;
      e_parity = (cd->ddom ? 1 : 0);
      curr_parity = ( e_parity == curr_parity ? 0 : 1 );
      rem_parity = ( st_parity ? 1 - curr_parity : curr_parity );

/*
      fprintf(stderr, "(%d,%d) [%d,%lf] ", curr_e->tail->id, curr_e->head->id,
              e_parity, EGdijkstraCostToLf(cd->weight));
*/

      e_val = cd->weight;
      curr_val = EGdijkstraCostAdd(curr_val, e_val);
      rem_ubound = EGdijkstraCostSub(st_ubound, curr_val);
 
      curr_s = curr_e->head;
      EGlistPushBack(curr_path, curr_e);
 
      /* check if we reached the end */
      if ( (curr_s->id/2) == (t->id/2) )
        break;
    }

    //fprintf(stderr, "\n");  

    /* check if we looped out */
    if ( cnt == gdata->cycleG->nodes->size )
    {
      EGlistClear(curr_path, nullFree);
      MESSAGE(1, "WARNING: GOT CAUGHT IN AN INFINITE LOOP");
    }

    /* check parity of the path obtained */ 
    if ( curr_path->size && (curr_parity != st_parity) )
    {
        //fprintf(stderr, "bad parity.\n");
        EGlistClear(curr_path, nullFree);
    }

    /* keep the best list */
    if ( curr_path->size && EGdijkstraCostIsLess(curr_val, *path_val) )
    {
      EGlistClear(path,nullFree);
      while(curr_path->size)
      { 
        EGlistPushBack(path, curr_path->begin->this);
        EGlistErase(curr_path, curr_path->begin, nullFree);
      }
      *path_val = curr_val;
    }
    EGlistClear(curr_path, nullFree);

  }

  EGfreeList(curr_path);
  
  return 0;

}

EGdGraphEdge_t* EGgreedySampleEdge(EGdGraphNode_t *s, 
                                   EGdGraphNode_t *t, 
                                   EGdijkstraCost_t st_ubound,
                                   int st_parity,
                                   EGgreedyData_t *gdata)
{
 
  int i, e_parity, cnt; 
  unsigned int bdg_head, bdg_t;
  EGlistNode_t *it;
  EGdGraphEdge_t *next_e=0, *e;
  EGdijkstraCost_t slack, st_lbound;
  double tot = 0.0;
  double rvalue;
  EGcycleData_t *cd;

  double *cand_weight = gdata->pedgevals;
  EGdGraphEdge_t **cand_edge = gdata->pedges;

  //fprintf(stderr, "st_ubound = %lf, st_parity = %d\n", EGdijkstraCostToLf(st_ubound), st_parity);

  //fprintf(stderr," s_out_edges = %u\n", s->out_edges->size);

  for(i=0, cnt=0, it=s->out_edges->begin; it && cnt < gdata->norig_edges; it=it->next, i++)
  {

    e = (EGdGraphEdge_t*)(it->this);
    cd = (EGcycleData_t*)(e->data);

    e_parity = (cd->ddom ? 1 : 0);

    bdg_head = gdata->bdGnodeMap[ e->head->id / 2 ]->id;
    bdg_t = gdata->bdGnodeMap[ t->id / 2 ]->id;

    // check if we have visitied the head node already
    if ( gdata->bdnodes[bdg_head] )
      continue;

    // we add an extra comparison because even/odd distance might == EG_DIJKSTRA_COST_MAX
    if ( st_parity == e_parity )
    {
      if (EGdijkstraCostIsLess(st_ubound, EGgetEvenDistance(bdg_head, bdg_t, gdata)))
        continue;
      st_lbound = EGdijkstraCostAdd( cd->weight, EGgetEvenDistance(bdg_head, bdg_t, gdata) );
    }
    else
    {
      if (EGdijkstraCostIsLess(st_ubound, EGgetOddDistance(bdg_head, bdg_t, gdata)))
        continue;
      st_lbound = EGdijkstraCostAdd( cd->weight, EGgetOddDistance(bdg_head, bdg_t, gdata) );
    }

    slack = EGdijkstraCostSub( st_ubound, st_lbound );

    #if 0
    fprintf(stderr, "%d candidate: (%d, %d) [%d, %lf] dist(%d, %d, %d) = %lf slack=%lf\n", i, 
                     e->tail->id, e->head->id, e_parity, EGdijkstraCostToLf(cd->weight),
                     e->tail->id, e->head->id, t->id, EGdijkstraCostToLf(st_lbound), EGdijkstraCostToLf(slack));
    #endif

    if ( EGdijkstraCostIsLess(slack,EGdijkstraToCost(0.0)) )
      continue;

    cand_edge[cnt] = e;
    cand_weight[cnt] = EGdijkstraCostToLf(slack);
    tot += cand_weight[cnt];
    cnt += 1;

  } 

  if (!cnt || cnt == gdata->norig_edges)
    goto DONE;

  next_e = cand_edge[0];

  rvalue = random();
  rvalue /= EGRAND_MAX;
  rvalue *= tot;

  tot = 0.0;

  for(i=0; i<cnt; i++)
  {
    tot += cand_weight[i];
    if ( rvalue < tot )
    {
      next_e = cand_edge[i];
      break;
    }
  }

  DONE:

  return next_e;

}

// assumes the list is sorted, from least slack, to most slack
// clears the item if it is not to go in
int EGaddIPtoList(EGinternalPairs_t *p, EGlist_t *pairs, unsigned int max_size)
{

  EGinternalPairs_t *q, *r;
  EGlistNode_t *it=0, *it2=0;

  /* if the list pairs is full, check to see if we have more slack 
     than every member of the list */
  if ( pairs->size && pairs->size == max_size )
  {
    q = (EGinternalPairs_t*)(pairs->end->this);
    if ( !EGdijkstraCostIsLess(p->value, q->value) )
      goto NO_ADD;
  }

  /* find the first element in the list who is not better than p */
  for(it = pairs->begin; it; it=it->next)
  {
    q = (EGinternalPairs_t*)(it->this); 
    if ( !EGdijkstraCostIsLess(q->value, p->value) )
      break;
  }
  //TEST(!it, "bad list");

  /* make sure p is not already in the list */
  /* for this, check every r having the same slack as q */
  for(it2 = it; it2; it2=it2->next) 
  {
    r = (EGinternalPairs_t*)(it2->this);
    if (r->value == p->value)
    {
      if (r->s == p->s)
        if (r->t == p->t)
          if (r->extpath->size == p->extpath->size)
            if (r->extpath->begin->this == p->extpath->begin->this)
              if (r->extpath->end->this == p->extpath->end->this)
                goto NO_ADD;
    }
    else
     break;
  }

  /* insert p before q */
  if (it)
    EGlistInsertBefore(pairs,it,p);
  else
    EGlistPushBack(pairs, p);

  while(pairs->size > max_size)
    EGlistEraseMP(pairs, pairs->end, EGfreeInternalPairs, pairs->mempool);

  return 0;

  NO_ADD:

  EGfreeInternalPairs(p,pairs->mempool);

  return 0;

}

// careful! This frees the two incoming lists!!
EGlist_t* EGmergeIPlists(EGlist_t *first, EGlist_t *second)
{

  EGlist_t *mlist;
  EGlistNode_t *it1, *it2;
  EGinternalPairs_t *p1, *p2;

  mlist = EGnewList(first->mempool);

  for(it1=first->begin, it2=second->begin; it1 || it2; )
  {
    p1 = it1 ? (EGinternalPairs_t*)(it1->this) : 0; 
    p2 = it2 ? (EGinternalPairs_t*)(it2->this) : 0;
    if (!it1)
    {
      EGlistPushBack(mlist, p2);
      it2=it2->next;
    }
    else if (!it2)
    {
      EGlistPushBack(mlist, p1);
      it1=it1->next;
    }
    else if (p1->value < p2->value)
    {
      EGlistPushBack(mlist, p1);
      it1=it1->next;
    }
    else
    {
      EGlistPushBack(mlist, p2);
      it2=it2->next;
    }
  }

  EGlistClear(first, nullFree);
  EGlistClear(second, nullFree);
  EGfreeList(first);
  EGfreeList(second);

  return mlist;

}

---