dbl_price.c

Go to the documentation of this file.

/****************************************************************************/
/*                                                                          */
/*  This file is part of QSopt_ex.                                          */
/*                                                                          */
/*  (c) Copyright 2006 by David Applegate, William Cook, Sanjeeb Dash,      */
/*  and Daniel Espinoza                                                     */
/*                                                                          */
/*  Sanjeeb Dash ownership of copyright in QSopt_ex is derived from his     */
/*  copyright in QSopt.                                                     */
/*                                                                          */
/*  This code may be used under the terms of the GNU General Public License */
/*  (Version 2.1 or later) as published by the Free Software Foundation.    */
/*                                                                          */
/*  Alternatively, use is granted for research purposes only.               */
/*                                                                          */
/*  It is your choice of which of these two licenses you are operating      */
/*  under.                                                                  */
/*                                                                          */
/*  We make no guarantees about the correctness or usefulness of this code. */
/*                                                                          */
/****************************************************************************/

/* RCS_INFO = "$RCSfile: dbl_price.c,v $ $Revision: 1.2 $ $Date: 2003/11/05 16:49:52 $"; */
//static int TRACE = 0;
#include "econfig.h"
#include "stddefs.h"
#include "dbl_qsopt.h"
#include "dbl_lpdefs.h"
#include "dbl_fct.h"
#include "dbl_price.h"
#include "dbl_basis.h"
#include "dbl_iqsutil.h"
#include "dbl_dstruct.h"
#ifdef USEDMALLOC
#include "dmalloc.h"
#endif

#define  dbl_MULTIP 1
#define  dbl_PRICE_DEBUG 0

static void dbl_update_d_scaleinf (
  dbl_price_info * const p,
  dbl_heap * const h,
  int const j,
  double inf,
  int const prule),
  dbl_update_p_scaleinf (
  dbl_price_info * const p,
  dbl_heap * const h,
  int const i,
  double inf,
  int const prule);

static void dbl_compute_dualI_inf (
  dbl_lpinfo * const lp,
  int const j,
  double * const inf),
  dbl_compute_dualII_inf (
  dbl_lpinfo * const lp,
  int const j,
  double * const inf),
  dbl_compute_primalI_inf (
  dbl_lpinfo * const lp,
  int const i,
  double * const inf),
  dbl_compute_primalII_inf (
  dbl_lpinfo * const lp,
  int const i,
  double * const inf);

void dbl_ILLprice_free_heap (
  dbl_price_info * const pinf)
{
  dbl_ILLheap_free (&(pinf->h));
}

int dbl_ILLprice_build_heap (
  dbl_price_info * const pinf,
  int const nkeys,
  double * keylist)
{
  dbl_ILLheap_init (&(pinf->h));
  dbl_EGlpNumSet (pinf->htrigger,
              1.0 +
              (double) nkeys / (PARAM_HEAP_RATIO * dbl_ILLutil_our_log2 (nkeys)));
  return dbl_ILLheap_build (&(pinf->h), nkeys, keylist);
}

int dbl_ILLprice_test_for_heap (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const nkeys,
  double * keylist,
  int const algo,
  int const upd)
{
  dbl_heap *const h = &(pinf->h);
  int rval = 0;
  double ravg;

  if (upd != 0)
  {
    dbl_EGlpNumInitVar (ravg);
    if (algo == PRIMAL_SIMPLEX)
      dbl_EGlpNumCopy (ravg, lp->cnts->za_ravg);
    else
      dbl_EGlpNumCopy (ravg, lp->cnts->y_ravg);
    if (dbl_EGlpNumIsLeq (ravg, pinf->htrigger))
      pinf->hineff--;
    else
    {
      dbl_EGlpNumDivUiTo (ravg, 2U);
      if (dbl_EGlpNumIsLess (pinf->htrigger, ravg))
        pinf->hineff++;
    }
    dbl_EGlpNumClearVar (ravg);
  }
  if (h->hexist == 0 && pinf->hineff <= 0)
  {
    rval = dbl_ILLprice_build_heap (pinf, nkeys, keylist);
    CHECKRVALG (rval, CLEANUP);
  }
  else if (h->hexist != 0 && pinf->hineff >= PARAM_HEAP_UTRIGGER)
  {
    dbl_ILLprice_free_heap (pinf);
    /*
     * printf ("freeing dbl_heap ..\n");
     * printf ("iter = %d, ravg = %.2f, trigger = %.2f\n",
     * lp->cnts->tot_iter, ravg, pinf->htrigger);
     */
  }

CLEANUP:
  if (rval)
    dbl_ILLprice_free_heap (pinf);
  return rval;
}

void dbl_ILLprice_init_pricing_info (
  dbl_price_info * const pinf)
{
  pinf->p_strategy = -1;
  pinf->d_strategy = -1;
  pinf->pI_price = -1;
  pinf->pII_price = -1;
  pinf->dI_price = -1;
  pinf->dII_price = -1;
  pinf->cur_price = -1;
  pinf->p_scaleinf = 0;
  pinf->d_scaleinf = 0;
  pinf->pdinfo.norms = 0;
  pinf->pdinfo.refframe = 0;
  pinf->psinfo.norms = 0;
  pinf->ddinfo.norms = 0;
  pinf->ddinfo.refframe = 0;
  pinf->dsinfo.norms = 0;
  pinf->dmpinfo.gstart = pinf->pmpinfo.gstart = 0;
  pinf->dmpinfo.gshift = pinf->pmpinfo.gshift = 0;
  pinf->dmpinfo.gsize = pinf->pmpinfo.gsize = 0;
  pinf->dmpinfo.bucket = pinf->pmpinfo.bucket = 0;
  pinf->dmpinfo.perm = pinf->pmpinfo.perm = 0;
  pinf->dmpinfo.infeas = pinf->pmpinfo.infeas = 0;
  dbl_ILLheap_init (&(pinf->h));
  dbl_EGlpNumZero (pinf->htrigger);
  pinf->hineff = 0;
}

void dbl_ILLprice_free_pricing_info (
  dbl_price_info * const pinf)
{
  dbl_EGlpNumFreeArray (pinf->p_scaleinf);
  dbl_EGlpNumFreeArray (pinf->d_scaleinf);
  dbl_EGlpNumFreeArray (pinf->pdinfo.norms);
  ILL_IFFREE (pinf->pdinfo.refframe, int);
  dbl_EGlpNumFreeArray (pinf->psinfo.norms);
  dbl_EGlpNumFreeArray (pinf->ddinfo.norms);
  ILL_IFFREE (pinf->ddinfo.refframe, int);
  dbl_EGlpNumFreeArray (pinf->dsinfo.norms);

  dbl_ILLprice_free_mpartial_info (&(pinf->pmpinfo));
  dbl_ILLprice_free_mpartial_info (&(pinf->dmpinfo));
  dbl_ILLprice_free_heap (pinf);
}

int dbl_ILLprice_build_pricing_info (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const phase)
{
  int rval = 0;
  int p_price = -1;
  int d_price = -1;

  switch (phase)
  {
  case PRIMAL_PHASEI:
    p_price = pinf->pI_price;
    break;
  case PRIMAL_PHASEII:
    p_price = pinf->pII_price;
    break;
  case DUAL_PHASEI:
    d_price = pinf->dI_price;
    break;
  case DUAL_PHASEII:
    d_price = pinf->dII_price;
    break;
  }

  if (p_price != -1)
  {
    pinf->cur_price = p_price;

    if (p_price == QS_PRICE_PDANTZIG || p_price == QS_PRICE_PDEVEX ||
        p_price == QS_PRICE_PSTEEP)
    {
      pinf->p_strategy = COMPLETE_PRICING;
      dbl_EGlpNumFreeArray (pinf->d_scaleinf);
      pinf->d_scaleinf = dbl_EGlpNumAllocArray (lp->nnbasic);
    }
    else if (p_price == QS_PRICE_PMULTPARTIAL)
      pinf->p_strategy = MULTI_PART_PRICING;

    switch (p_price)
    {
    case QS_PRICE_PDEVEX:
      if (pinf->pdinfo.norms)
        return rval;
      rval = dbl_ILLprice_build_pdevex_norms (lp, &(pinf->pdinfo), 0);
      CHECKRVALG(rval,CLEANUP);
      break;
    case QS_PRICE_PSTEEP:
      if (pinf->psinfo.norms)
        return rval;
      rval = dbl_ILLprice_build_psteep_norms (lp, &(pinf->psinfo));
      CHECKRVALG(rval,CLEANUP);
      break;
    case QS_PRICE_PMULTPARTIAL:
      rval = dbl_ILLprice_build_mpartial_info (lp, pinf, COL_PRICING);
      CHECKRVALG(rval,CLEANUP);
      break;
    }
  }
  else if (d_price != -1)
  {
    pinf->cur_price = d_price;

    if (d_price == QS_PRICE_DDANTZIG || d_price == QS_PRICE_DSTEEP ||
        d_price == QS_PRICE_DDEVEX)
    {
      pinf->d_strategy = COMPLETE_PRICING;
      dbl_EGlpNumFreeArray (pinf->p_scaleinf);
      pinf->p_scaleinf = dbl_EGlpNumAllocArray (lp->nrows);
    }
    else if (d_price == QS_PRICE_DMULTPARTIAL)
      pinf->d_strategy = MULTI_PART_PRICING;

    switch (d_price)
    {
    case QS_PRICE_DSTEEP:
      if (pinf->dsinfo.norms)
        return rval;
      rval = dbl_ILLprice_build_dsteep_norms (lp, &(pinf->dsinfo));
      CHECKRVALG(rval,CLEANUP);
      break;
    case QS_PRICE_DMULTPARTIAL:
      rval = dbl_ILLprice_build_mpartial_info (lp, pinf, ROW_PRICING);
      CHECKRVALG(rval,CLEANUP);
      break;
    case QS_PRICE_DDEVEX:
      if (pinf->ddinfo.norms)
        return rval;
      rval = dbl_ILLprice_build_ddevex_norms (lp, &(pinf->ddinfo), 0);
      CHECKRVALG(rval,CLEANUP);
      break;
    }
  }

CLEANUP:
  if (rval)
    dbl_ILLprice_free_pricing_info (pinf);
  EG_RETURN(rval);
}

int dbl_ILLprice_update_pricing_info (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const phase,
  dbl_svector * const wz,
  int const eindex,
  int const lindex,
  double y)
{
  int rval = 0;
  int p_price = -1;
  int d_price = -1;

  switch (phase)
  {
  case PRIMAL_PHASEI:
    p_price = pinf->pI_price;
    break;
  case PRIMAL_PHASEII:
    p_price = pinf->pII_price;
    break;
  case DUAL_PHASEI:
    d_price = pinf->dI_price;
    break;
  case DUAL_PHASEII:
    d_price = pinf->dII_price;
    break;
  }

  if (p_price != -1)
  {
    if (p_price == QS_PRICE_PDEVEX)
    {
      rval = dbl_ILLprice_update_pdevex_norms (lp, &(pinf->pdinfo), eindex, y);
      CHECKRVALG(rval,CLEANUP);
    }
    else if (p_price == QS_PRICE_PSTEEP)
      dbl_ILLprice_update_psteep_norms (lp, &(pinf->psinfo), wz, eindex, y);
  }
  else if (d_price != -1)
  {
    if (d_price == QS_PRICE_DSTEEP)
      dbl_ILLprice_update_dsteep_norms (lp, &(pinf->dsinfo), wz, lindex, y);
    else if (d_price == QS_PRICE_DDEVEX)
    {
      rval = dbl_ILLprice_update_ddevex_norms (lp, &(pinf->ddinfo), lindex, y);
      CHECKRVALG(rval,CLEANUP);
    }
  }
CLEANUP:
  EG_RETURN(rval);
}

int dbl_ILLprice_get_price (
  dbl_price_info * const p,
  int const phase)
{
  int pri = -1;

  switch (phase)
  {
  case PRIMAL_PHASEI:
    return p->pI_price;
  case PRIMAL_PHASEII:
    return p->pII_price;
  case DUAL_PHASEI:
    return p->dI_price;
  case DUAL_PHASEII:
    return p->dII_price;
  }
  return pri;
}

void dbl_ILLprice_free_mpartial_info (
  dbl_mpart_info * p)
{
  ILL_IFFREE (p->gstart, int);
  ILL_IFFREE (p->gshift, int);
  ILL_IFFREE (p->gsize, int);
  ILL_IFFREE (p->bucket, int);
  dbl_EGlpNumFreeArray (p->infeas);
  ILL_IFFREE (p->perm, int);
}

int dbl_ILLprice_build_mpartial_info (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const pricetype)
{
  int i = 0;
  int rval = 0;
  int extra = 0;
  int nelems;
  dbl_mpart_info *p;

  p = (pricetype == COL_PRICING) ? &(pinf->pmpinfo) : &(pinf->dmpinfo);
  p->k = 50;
  p->cgroup = 0;
  nelems = (pricetype == COL_PRICING) ? lp->nnbasic : lp->nrows;

  if (nelems % p->k)
    extra = nelems - p->k * (nelems / p->k);
  p->ngroups = nelems / p->k;
  if (extra != 0)
    p->ngroups++;

  ILL_SAFE_MALLOC (p->gstart, p->ngroups, int);
  ILL_SAFE_MALLOC (p->gshift, p->ngroups, int);
  ILL_SAFE_MALLOC (p->gsize, p->ngroups, int);
  ILL_SAFE_MALLOC (p->bucket, 2 * p->k, int);
  p->infeas = dbl_EGlpNumAllocArray (2 * p->k);
  ILL_SAFE_MALLOC (p->perm, 2 * p->k, int);

  p->bsize = 0;

  if (extra != 0)
  {
    p->gstart[0] = 0;
    p->gshift[0] = 1;
    p->gsize[0] = extra;
    for (i = 1; i < p->ngroups; i++)
    {
      p->gstart[i] = extra + i - 1;
      p->gshift[i] = p->ngroups - 1;
      p->gsize[i] = p->k;
    }
  }
  else
  {
    for (i = 0; i < p->ngroups; i++)
    {
      p->gstart[i] = i;
      p->gshift[i] = p->ngroups;
      p->gsize[i] = p->k;
    }
  }

CLEANUP:
  if (rval)
    dbl_ILLprice_free_mpartial_info (p);
  EG_RETURN(rval);
}

void dbl_ILLprice_init_mpartial_price (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const phase,
  int const pricetype)
{
  int i;
  dbl_mpart_info *p;

  p = (pricetype == COL_PRICING) ? &(pinf->pmpinfo) : &(pinf->dmpinfo);
  p->bsize = 0;
  i = p->cgroup;
  do
  {
    dbl_ILLprice_mpartial_group (lp, p, phase, i, pricetype);
    i = (i + 1) % p->ngroups;
  } while (i != p->cgroup && p->bsize <= p->k);
  p->cgroup = i;
}

void dbl_ILLprice_update_mpartial_price (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const phase,
  int const pricetype)
{
  int i = 0;
  int csize = 0;
  double infeas;
  dbl_mpart_info *p;
  dbl_price_res pr;

  p = (pricetype == COL_PRICING) ? &(pinf->pmpinfo) : &(pinf->dmpinfo);
  dbl_EGlpNumInitVar (pr.dinfeas);
  dbl_EGlpNumInitVar (pr.pinfeas);
  dbl_EGlpNumInitVar (infeas);

#ifdef dbl_MULTIP
  i = 0;
  while (i < p->bsize)
  {
    if (pricetype == COL_PRICING)
    {
      dbl_ILLprice_column (lp, p->bucket[i], phase, &pr);
      dbl_EGlpNumCopy (infeas, pr.dinfeas);
    }
    else
    {
      dbl_ILLprice_row (lp, p->bucket[i], phase, &pr);
      dbl_EGlpNumCopy (infeas, pr.pinfeas);
    }
    if (!dbl_EGlpNumIsNeqqZero (infeas))
    {
      p->bucket[i] = p->bucket[p->bsize - 1];
      p->bsize--;
    }
    else
    {
      dbl_EGlpNumCopy (p->infeas[i], infeas);
      i++;
    }
  }
  if (p->bsize > 0)
  {
    for (i = 0; i < p->bsize; i++)
      p->perm[i] = i;
    dbl_EGutilPermSort ((size_t) (p->bsize), p->perm,
                    (const double * const) p->infeas);

    csize = MIN (p->bsize, p->k);
    for (i = csize - 1; i >= 0; i--)
      lp->iwork[p->bucket[p->perm[i]]] = 1;

    for (i = 0, csize = 0; i < p->bsize; i++)
      if (lp->iwork[p->bucket[i]] == 1)
      {
        dbl_EGlpNumCopy (p->infeas[csize], p->infeas[i]);
        p->bucket[csize] = p->bucket[i];
        csize++;
      }
    p->bsize = csize;
  }
#else
  p->bsize = 0;
#endif

  i = p->cgroup;
  do
  {
    dbl_ILLprice_mpartial_group (lp, p, phase, i, pricetype);
    i = (i + 1) % p->ngroups;
  } while (i != p->cgroup && p->bsize <= p->k);
  p->cgroup = i;

#ifdef dbl_MULTIP
  for (i = 0; i < csize; i++)
    lp->iwork[p->bucket[i]] = 0;
#endif
  dbl_EGlpNumClearVar (infeas);
  dbl_EGlpNumClearVar (pr.pinfeas);
  dbl_EGlpNumClearVar (pr.dinfeas);
}

void dbl_ILLprice_delete_onempart_price (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const indx,
  int const pricetype)
{
  int i = 0;
  dbl_mpart_info *p;

  p = (pricetype == COL_PRICING) ? &(pinf->pmpinfo) : &(pinf->dmpinfo);

  for (i = 0; i < p->bsize; i++)
    if (p->bucket[i] == indx)
    {
      p->bucket[i] = p->bucket[p->bsize - 1];
      dbl_EGlpNumCopy (p->infeas[i], p->infeas[p->bsize - 1]);
      p->bsize--;
      break;
    }
}

void dbl_ILLprice_mpartial_group (
  dbl_lpinfo * const lp,
  dbl_mpart_info * const p,
  int const phase,
  int const g,
  int const pricetype)
{
  int i, ix;
  int gstart = p->gstart[g];
  int gsize = p->gsize[g];
  int gshift = p->gshift[g];
  double infeas;
  dbl_price_res pr;

  dbl_EGlpNumInitVar (pr.dinfeas);
  dbl_EGlpNumInitVar (pr.pinfeas);
  dbl_EGlpNumInitVar (infeas);

  for (i = 0, ix = gstart; i < gsize; i++, ix += gshift)
  {
#ifdef dbl_MULTIP
    if (lp->iwork[ix])
      continue;
#endif
    if (pricetype == COL_PRICING)
    {
      dbl_ILLprice_column (lp, ix, phase, &pr);
      dbl_EGlpNumCopy (infeas, pr.dinfeas);
    }
    else
    {
      dbl_ILLprice_row (lp, ix, phase, &pr);
      dbl_EGlpNumCopy (infeas, pr.pinfeas);
    }
    if (dbl_EGlpNumIsNeqqZero (infeas))
    {
      dbl_EGlpNumCopy (p->infeas[p->bsize], infeas);
      p->bucket[p->bsize] = ix;
      p->bsize++;
    }
  }
  dbl_EGlpNumClearVar (infeas);
  dbl_EGlpNumClearVar (pr.dinfeas);
  dbl_EGlpNumClearVar (pr.pinfeas);
}

void dbl_ILLprice_column (
  dbl_lpinfo * const lp,
  int const ix,
  int const phase,
  dbl_price_res * const pr)
{
  int i;
  int col;
  int mcnt;
  int mbeg;
  double sum;

  dbl_EGlpNumZero (pr->dinfeas);
  col = lp->nbaz[ix];
  if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED)
    return;
  dbl_EGlpNumInitVar (sum);
  dbl_EGlpNumZero (sum);
  mcnt = lp->matcnt[col];
  mbeg = lp->matbeg[col];

  if (phase == PRIMAL_PHASEII)
  {
    for (i = 0; i < mcnt; i++)
      dbl_EGlpNumAddInnProdTo (sum, lp->piz[lp->matind[mbeg + i]], lp->matval[mbeg + i]);
    dbl_EGlpNumCopyDiff (lp->dz[ix], lp->cz[col], sum);
    dbl_compute_dualII_inf (lp, ix, &(pr->dinfeas));
  }
  else
  {
    for (i = 0; i < mcnt; i++)
      dbl_EGlpNumAddInnProdTo (sum, lp->pIpiz[lp->matind[mbeg + i]], lp->matval[mbeg + i]);
    dbl_EGlpNumCopyNeg (lp->pIdz[ix], sum);
    dbl_compute_dualI_inf (lp, ix, &(pr->dinfeas));
  }
  dbl_EGlpNumClearVar (sum);
}

void dbl_ILLprice_row (
  dbl_lpinfo * const lp,
  int const ix,
  int const phase,
  dbl_price_res * const pr)
{
  if (phase == DUAL_PHASEII)
    dbl_compute_primalII_inf (lp, ix, &(pr->pinfeas));
  else
    dbl_compute_primalI_inf (lp, ix, &(pr->pinfeas));
}

int dbl_ILLprice_build_pdevex_norms (
  dbl_lpinfo * const lp,
  dbl_p_devex_info * const pdinfo,
  int const reinit)
{
  int j;
  int rval = 0;

  if (reinit == 0)
  {
    pdinfo->ninit = 0;
    pdinfo->norms = dbl_EGlpNumAllocArray (lp->nnbasic);
    ILL_SAFE_MALLOC (pdinfo->refframe, lp->ncols, int);
  }

  if (reinit != 0)
    pdinfo->ninit++;

  for (j = 0; j < lp->ncols; j++)
  {
    if (lp->vstat[j] == STAT_BASIC)
      pdinfo->refframe[j] = 0;
    else
    {
      dbl_EGlpNumOne (pdinfo->norms[lp->vindex[j]]);
      pdinfo->refframe[j] = 1;
    }
  }

CLEANUP:
  if (rval)
  {
    dbl_EGlpNumFreeArray (pdinfo->norms);
    ILL_IFFREE (pdinfo->refframe, int);
  }
  EG_RETURN(rval);
}

int dbl_ILLprice_update_pdevex_norms (
  dbl_lpinfo * const lp,
  dbl_p_devex_info * const pdinfo,
  int const eindex,
  double yl)
{
  int i, j;
  double normj;
  double zAj;
  double ntmp, ntmp2;

  dbl_EGlpNumInitVar (normj);
  dbl_EGlpNumInitVar (zAj);
  dbl_EGlpNumInitVar (ntmp);
  dbl_EGlpNumInitVar (ntmp2);
  dbl_EGlpNumZero (normj);

  for (i = 0; i < lp->yjz.nzcnt; i++)
    if (pdinfo->refframe[lp->baz[lp->yjz.indx[i]]])
      dbl_EGlpNumAddInnProdTo (normj, lp->yjz.coef[i], lp->yjz.coef[i]);

  if (pdinfo->refframe[lp->nbaz[eindex]])
    dbl_EGlpNumAddTo (normj, dbl_oneLpNum);

  dbl_EGlpNumSet(ntmp,1000.0);
  dbl_EGlpNumSet(ntmp2,0.001);
  dbl_EGlpNumMultTo(ntmp,pdinfo->norms[eindex]);
  dbl_EGlpNumMultTo(ntmp2,pdinfo->norms[eindex]);
  if (dbl_EGlpNumIsLess (normj, ntmp2) || dbl_EGlpNumIsLess (ntmp, normj))
  {
    dbl_EGlpNumClearVar (zAj);
    dbl_EGlpNumClearVar (normj);
    dbl_EGlpNumClearVar (ntmp);
    dbl_EGlpNumClearVar (ntmp2);
    return dbl_ILLprice_build_pdevex_norms (lp, pdinfo, 1);
  }

  for (i = 0; i < lp->zA.nzcnt; i++)
  {
    j = lp->zA.indx[i];
    dbl_EGlpNumCopyFrac (zAj, lp->zA.coef[i], yl);
    dbl_EGlpNumMultTo (zAj, zAj);
    dbl_EGlpNumMultTo (zAj, normj);
    if (dbl_EGlpNumIsLess (pdinfo->norms[j], zAj))
      dbl_EGlpNumCopy (pdinfo->norms[j], zAj);
  }
  dbl_EGlpNumDivTo (normj, yl);
  dbl_EGlpNumDivTo (normj, yl);
  if (dbl_EGlpNumIsLess (normj, dbl_oneLpNum))
    dbl_EGlpNumCopy (pdinfo->norms[eindex], dbl_oneLpNum);
  else
    dbl_EGlpNumCopy (pdinfo->norms[eindex], normj);
  dbl_EGlpNumClearVar (zAj);
  dbl_EGlpNumClearVar (normj);
  dbl_EGlpNumClearVar (ntmp);
  dbl_EGlpNumClearVar (ntmp2);
  return 0;
}

int dbl_ILLprice_build_psteep_norms (
  dbl_lpinfo * const lp,
  dbl_p_steep_info * const psinfo)
{
  int j;
  int rval = 0;
  dbl_svector yz;

  dbl_ILLsvector_init (&yz);
  rval = dbl_ILLsvector_alloc (&yz, lp->nrows);
  CHECKRVALG(rval,CLEANUP);
  psinfo->norms = dbl_EGlpNumAllocArray (lp->nnbasic);

  for (j = 0; j < lp->nnbasic; j++)
  {
    rval = ILLstring_report (NULL, &lp->O->reporter);
    CHECKRVALG(rval,CLEANUP);
    dbl_ILLfct_compute_yz (lp, &yz, 0, lp->nbaz[j]);
    dbl_EGlpNumInnProd (psinfo->norms[j], yz.coef, yz.coef, (size_t) yz.nzcnt);
    dbl_EGlpNumAddTo (psinfo->norms[j], dbl_oneLpNum);
  }

CLEANUP:
  dbl_ILLsvector_free (&yz);
  if (rval)
    dbl_EGlpNumFreeArray (psinfo->norms);

  EG_RETURN(rval);
}

void dbl_ILLprice_update_psteep_norms (
  dbl_lpinfo * const lp,
  dbl_p_steep_info * const psinfo,
  dbl_svector * const wz,
  int const eindex,
  double yl)
{
  int i, j, k;
  int mcnt, mbeg;
  double normj,ntmp;
  double zAj, wAj;
  double *v = 0;

  dbl_EGlpNumInitVar (normj);
  dbl_EGlpNumInitVar (zAj);
  dbl_EGlpNumInitVar (ntmp);
  dbl_EGlpNumInitVar (wAj);
  dbl_EGlpNumInnProd (normj, lp->yjz.coef, lp->yjz.coef, (size_t) (lp->yjz.nzcnt));
  dbl_EGlpNumAddTo (normj, dbl_oneLpNum);

#if 0
  Bico - remove warnings for dist
    if (fabs ((normj - psinfo->norms[eindex]) / normj) > 1000.0 /* 0.01 */ )
    {
      printf ("warning: incorrect norm values\n");
      printf ("anorm = %.6f, pnorm = %.6f\n", normj, psinfo->norms[eindex]);
      fflush (stdout);
    }
#endif

  dbl_ILLfct_load_workvector (lp, wz);
  v = lp->work.coef;

  for (k = 0; k < lp->zA.nzcnt; k++)
  {
    j = lp->zA.indx[k];
    dbl_EGlpNumCopy (zAj, lp->zA.coef[k]);
    dbl_EGlpNumZero (wAj);
    mcnt = lp->matcnt[lp->nbaz[j]];
    mbeg = lp->matbeg[lp->nbaz[j]];
    for (i = 0; i < mcnt; i++)
      dbl_EGlpNumAddInnProdTo (wAj, lp->matval[mbeg + i], v[lp->matind[mbeg + i]]);

    /* compute ntmp = (zAj * ((zAj * normj / yl) - (2.0 * wAj))) / yl; */ 
    dbl_EGlpNumCopy(ntmp,zAj);
    dbl_EGlpNumMultTo(ntmp,normj);
    dbl_EGlpNumDivTo(ntmp,yl);
    dbl_EGlpNumSubTo(ntmp,wAj);
    dbl_EGlpNumSubTo(ntmp,wAj);
    dbl_EGlpNumMultTo(ntmp,zAj);
    dbl_EGlpNumDivTo(ntmp,yl);
    /* set psinfo->norms[j] += (zAj * ((zAj * normj / yl) - (2.0 * wAj))) / yl; */
    dbl_EGlpNumAddTo(psinfo->norms[j],ntmp);
    if (dbl_EGlpNumIsLess (psinfo->norms[j], dbl_oneLpNum))
      dbl_EGlpNumOne (psinfo->norms[j]);
  }

  dbl_EGlpNumCopyFrac (psinfo->norms[eindex], normj, yl);
  dbl_EGlpNumDivTo (psinfo->norms[eindex], yl);
  if (dbl_EGlpNumIsLess (psinfo->norms[eindex], dbl_oneLpNum))
    dbl_EGlpNumOne (psinfo->norms[eindex]);

  dbl_ILLfct_zero_workvector (lp);
  dbl_EGlpNumClearVar (wAj);
  dbl_EGlpNumClearVar (zAj);
  dbl_EGlpNumClearVar (normj);
  dbl_EGlpNumClearVar (ntmp);
}

int dbl_ILLprice_build_ddevex_norms (
  dbl_lpinfo * const lp,
  dbl_d_devex_info * const ddinfo,
  int const reinit)
{
  int i;
  int rval = 0;

  if (reinit == 0)
  {
    ddinfo->ninit = 0;
    ddinfo->norms = dbl_EGlpNumAllocArray (lp->nrows);
    ILL_SAFE_MALLOC (ddinfo->refframe, lp->ncols, int);
  }
  if (reinit != 0)
    ddinfo->ninit++;

  for (i = 0; i < lp->ncols; i++)
    ddinfo->refframe[i] = (lp->vstat[i] == STAT_BASIC) ? 1 : 0;

  for (i = 0; i < lp->nrows; i++)
    dbl_EGlpNumOne (ddinfo->norms[i]);

CLEANUP:
  if (rval)
  {
    dbl_EGlpNumFreeArray (ddinfo->norms);
    ILL_IFFREE (ddinfo->refframe, int);
  }
  EG_RETURN(rval);
}

int dbl_ILLprice_update_ddevex_norms (
  dbl_lpinfo * const lp,
  dbl_d_devex_info * const ddinfo,
  int const lindex,
  double yl)
{
  int i, r;
  double normi;
  double yr;
  double ntmp,ntmp2;

  dbl_EGlpNumInitVar (ntmp);
  dbl_EGlpNumInitVar (ntmp2);
  dbl_EGlpNumInitVar (normi);
  dbl_EGlpNumInitVar (yr);
  dbl_EGlpNumZero (normi);

  for (i = 0; i < lp->zA.nzcnt; i++)
    if (ddinfo->refframe[lp->nbaz[lp->zA.indx[i]]])
      dbl_EGlpNumAddInnProdTo (normi, lp->zA.coef[i], lp->zA.coef[i]);

  if (ddinfo->refframe[lp->baz[lindex]])
    dbl_EGlpNumAddTo (normi, dbl_oneLpNum);

  dbl_EGlpNumSet(ntmp,1000.0);
  dbl_EGlpNumSet(ntmp2,0.001);
  dbl_EGlpNumMultTo(ntmp,ddinfo->norms[lindex]);
  dbl_EGlpNumMultTo(ntmp2,ddinfo->norms[lindex]);
  if (dbl_EGlpNumIsLess(normi, ntmp2) || dbl_EGlpNumIsLess(ntmp, normi))
  {
    dbl_EGlpNumClearVar (normi);
    dbl_EGlpNumClearVar (yr);
    dbl_EGlpNumClearVar (ntmp);
    dbl_EGlpNumClearVar (ntmp2);
    return dbl_ILLprice_build_ddevex_norms (lp, ddinfo, 1);
  }

  for (i = 0; i < lp->yjz.nzcnt; i++)
  {
    r = lp->yjz.indx[i];
    dbl_EGlpNumCopy(yr, lp->yjz.coef[i]);
    dbl_EGlpNumCopy(ntmp,yr);
    dbl_EGlpNumMultTo(ntmp,yr);
    dbl_EGlpNumMultTo(ntmp,normi);
    dbl_EGlpNumDivTo(ntmp,yl);
    dbl_EGlpNumDivTo(ntmp,yl);
    if (dbl_EGlpNumIsLess (ddinfo->norms[r], ntmp))
      dbl_EGlpNumCopy (ddinfo->norms[r], ntmp);
  }
  dbl_EGlpNumCopy (ddinfo->norms[lindex], normi);
  dbl_EGlpNumDivTo(ddinfo->norms[lindex], yl);
  dbl_EGlpNumDivTo(ddinfo->norms[lindex], yl);
  if (dbl_EGlpNumIsLess (ddinfo->norms[lindex], dbl_oneLpNum))
    dbl_EGlpNumOne (ddinfo->norms[lindex]);
  dbl_EGlpNumClearVar (normi);
  dbl_EGlpNumClearVar (yr);
  dbl_EGlpNumClearVar (ntmp);
  dbl_EGlpNumClearVar (ntmp2);
  return 0;
}

int dbl_ILLprice_build_dsteep_norms (
  dbl_lpinfo * const lp,
  dbl_d_steep_info * const dsinfo)
{
  int i;
  int rval = 0;
  dbl_svector z;

  dbl_ILLsvector_init (&z);
  rval = dbl_ILLsvector_alloc (&z, lp->nrows);
  CHECKRVALG(rval,CLEANUP);
  dsinfo->norms = dbl_EGlpNumAllocArray (lp->nrows);

  for (i = 0; i < lp->nrows; i++)
  {
    rval = ILLstring_report (NULL, &lp->O->reporter);
    CHECKRVALG(rval,CLEANUP);

    dbl_ILLfct_compute_zz (lp, &z, i);

    dbl_EGlpNumInnProd (dsinfo->norms[i], z.coef, z.coef, (size_t) z.nzcnt);
    if (dbl_EGlpNumIsLess (dsinfo->norms[i], dbl_PARAM_MIN_DNORM))
      dbl_EGlpNumCopy (dsinfo->norms[i], dbl_PARAM_MIN_DNORM);
  }

CLEANUP:
  dbl_ILLsvector_free (&z);
  if (rval)
    dbl_EGlpNumFreeArray (dsinfo->norms);

  EG_RETURN(rval);
}

int dbl_ILLprice_get_dsteep_norms (
  dbl_lpinfo * const lp,
  int const count,
  int *const rowind,
  double * const norms)
{
  int i;
  int rval = 0;
  dbl_svector z;

  dbl_ILLsvector_init (&z);
  rval = dbl_ILLsvector_alloc (&z, lp->nrows);
  CHECKRVALG(rval,CLEANUP);

  for (i = 0; i < count; i++)
  {
    dbl_ILLfct_compute_zz (lp, &z, rowind[i]);
    dbl_EGlpNumInnProd (norms[i], z.coef, z.coef, (size_t) z.nzcnt);
  }

CLEANUP:
  dbl_ILLsvector_free (&z);
  EG_RETURN(rval);
}

void dbl_ILLprice_update_dsteep_norms (
  dbl_lpinfo * const lp,
  dbl_d_steep_info * const dsinfo,
  dbl_svector * const wz,
  int const lindex,
  double yl)
{
  int i, k;
  double yij;
  double norml;
  double *v = 0;
  double ntmp;

  dbl_EGlpNumInitVar (ntmp);
  dbl_EGlpNumInitVar (norml);
  dbl_EGlpNumInitVar (yij);
  dbl_EGlpNumInnProd (norml, lp->zz.coef, lp->zz.coef, (size_t) (lp->zz.nzcnt));

#if 0
  Bico - remove warnings for dist
    if (fabs ((norml - dsinfo->norms[lindex]) / norml) > 1000.0 /*0.01 */ )
    {
      printf ("warning: incorrect dnorm values\n");
      printf ("anorm = %.6f, pnorm = %.6f\n", norml, dsinfo->norms[lindex]);
      fflush (stdout);
    }
#endif

  dbl_ILLfct_load_workvector (lp, wz);
  v = lp->work.coef;

  for (k = 0; k < lp->yjz.nzcnt; k++)
  {
    i = lp->yjz.indx[k];
    dbl_EGlpNumCopy (yij, lp->yjz.coef[k]);
    /* compute in ntmp (yij * ((yij * norml / yl) - (2.0 * v[i]))) / yl; */
    dbl_EGlpNumCopy(ntmp,yij);
    dbl_EGlpNumMultTo(ntmp,norml);
    dbl_EGlpNumDivTo(ntmp,yl);
    dbl_EGlpNumSubTo(ntmp,v[i]);
    dbl_EGlpNumSubTo(ntmp,v[i]);
    dbl_EGlpNumMultTo (ntmp, yij);
    dbl_EGlpNumDivTo (ntmp, yl);
    /* set dsinfo->norms[i] += (yij * ((yij * norml / yl) - (2.0 * v[i]))) / yl;*/
    dbl_EGlpNumAddTo(dsinfo->norms[i], ntmp);
    if (dbl_EGlpNumIsLess (dsinfo->norms[i], dbl_PARAM_MIN_DNORM))
      dbl_EGlpNumCopy (dsinfo->norms[i], dbl_PARAM_MIN_DNORM);
  }
  dbl_EGlpNumCopyFrac (dsinfo->norms[lindex], norml, yl);
  dbl_EGlpNumDivTo (dsinfo->norms[lindex], yl);
  if (dbl_EGlpNumIsLess (dsinfo->norms[lindex], dbl_PARAM_MIN_DNORM))
    dbl_EGlpNumCopy (dsinfo->norms[lindex], dbl_PARAM_MIN_DNORM);

  dbl_ILLfct_zero_workvector (lp);
  dbl_EGlpNumClearVar (norml);
  dbl_EGlpNumClearVar (ntmp);
  dbl_EGlpNumClearVar (yij);
}

static void dbl_update_d_scaleinf (
  dbl_price_info * const p,
  dbl_heap * const h,
  int const j,
  double inf,
  int const prule)
{
  if (!dbl_EGlpNumIsNeqqZero (inf))
  {
    dbl_EGlpNumZero (p->d_scaleinf[j]);
    if (h->hexist != 0 && h->loc[j] != -1)
      dbl_ILLheap_delete (h, j);
  }
  else
  {
    if (prule == QS_PRICE_PDANTZIG)
      dbl_EGlpNumCopy (p->d_scaleinf[j], inf);
    else if (prule == QS_PRICE_PDEVEX)
      dbl_EGlpNumCopySqrOver (p->d_scaleinf[j], inf, p->pdinfo.norms[j]);
    else if (prule == QS_PRICE_PSTEEP)
      dbl_EGlpNumCopySqrOver (p->d_scaleinf[j], inf, p->psinfo.norms[j]);

    if (h->hexist != 0)
    {
      if (h->loc[j] == -1)
        dbl_ILLheap_insert (h, j);
      else
        dbl_ILLheap_modify (h, j);
    }
  }
}

static void dbl_compute_dualI_inf (
  dbl_lpinfo * const lp,
  const int j,
  double * const inf)
{
  int col = lp->nbaz[j];
  int vt = lp->vtype[col];
  int vs = lp->vstat[col];
  double*dj = &(lp->pIdz[j]);
  double*ftol = &(lp->tol->id_tol);
  dbl_EGlpNumZero (*inf);
  if (vt != VARTIFICIAL && vt != VFIXED)
  {
    if( dbl_EGlpNumIsSumLess(*dj,*ftol,dbl_zeroLpNum) && (vs == STAT_LOWER || vs == STAT_ZERO))
      dbl_EGlpNumCopyNeg(*inf,*dj);
    else if (dbl_EGlpNumIsLess(*ftol, *dj) && (vs == STAT_UPPER || vs == STAT_ZERO))
      dbl_EGlpNumCopy (*inf, *dj);
  }
}

static void dbl_compute_dualII_inf (
  dbl_lpinfo * const lp,
  int const j,
  double * const inf)
{
  int col = lp->nbaz[j];
  int vt = lp->vtype[col];
  int vs = lp->vstat[col];
  double*dj = &(lp->dz[j]);
  double*ftol = &(lp->tol->dfeas_tol);
  dbl_EGlpNumZero (*inf);
  if (vt != VARTIFICIAL && vt != VFIXED)
  {
    if( dbl_EGlpNumIsSumLess(*dj,*ftol,dbl_zeroLpNum) && (vs == STAT_LOWER || vs == STAT_ZERO))
      dbl_EGlpNumCopyNeg(*inf,*dj);
    else if (dbl_EGlpNumIsLess(*ftol,*dj) && (vs == STAT_UPPER || vs == STAT_ZERO))
      dbl_EGlpNumCopy (*inf, *dj);
  }
}

void dbl_ILLprice_compute_dual_inf (
  dbl_lpinfo * const lp,
  dbl_price_info * const p,
  int *const ix,
  int const icnt,
  int const phase)
{
  int i;
  int price;
  double inf;
  dbl_heap *h = &(p->h);

  price = (phase == PRIMAL_PHASEI) ? p->pI_price : p->pII_price;
  dbl_EGlpNumInitVar (inf);
  dbl_EGlpNumZero (inf);

  if (phase == PRIMAL_PHASEI)
  {
    if (ix == NULL)
      for (i = 0; i < lp->nnbasic; i++)
      {
        dbl_compute_dualI_inf (lp, i, &(inf));
        dbl_update_d_scaleinf (p, h, i, inf, price);
      }
    else
      for (i = 0; i < icnt; i++)
      {
        dbl_compute_dualI_inf (lp, ix[i], &(inf));
        dbl_update_d_scaleinf (p, h, ix[i], inf, price);
      }
  }
  else if (phase == PRIMAL_PHASEII)
  {
    if (ix == NULL)
      for (i = 0; i < lp->nnbasic; i++)
      {
        dbl_compute_dualII_inf (lp, i, &inf);
        dbl_update_d_scaleinf (p, h, i, inf, price);
      }
    else
      for (i = 0; i < icnt; i++)
      {
        dbl_compute_dualII_inf (lp, ix[i], &inf);
        dbl_update_d_scaleinf (p, h, ix[i], inf, price);
      }
  }
  dbl_EGlpNumClearVar (inf);
}

void dbl_ILLprice_primal (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  dbl_price_res * const pr,
  int const phase)
{
  int j, vs;
  double d_e, d_max;
  double *ftol = &(lp->tol->dfeas_tol);
  dbl_heap *const h = &(pinf->h);

  dbl_EGlpNumInitVar (d_e);
  dbl_EGlpNumInitVar (d_max);
  pr->eindex = -1;
  dbl_EGlpNumZero(d_max);

#if USEHEAP > 0
  dbl_ILLprice_test_for_heap (lp, pinf, lp->nnbasic, pinf->d_scaleinf,
                          PRIMAL_SIMPLEX, 1);
#endif

  if (pinf->p_strategy == COMPLETE_PRICING)
  {
    if (h->hexist)
    {
      pr->eindex = dbl_ILLheap_findmin (h);
      if (pr->eindex != -1)
        dbl_ILLheap_delete (h, pr->eindex);
    }
    else
    {
      for (j = 0; j < lp->nnbasic; j++)
      {
        if (dbl_EGlpNumIsLess (d_max, pinf->d_scaleinf[j]))
        {
          dbl_EGlpNumCopy (d_max, pinf->d_scaleinf[j]);
          pr->eindex = j;
        }
      }
    }
  }
  else if (pinf->p_strategy == MULTI_PART_PRICING)
  {
    for (j = 0; j < pinf->pmpinfo.bsize; j++)
    {
      if (dbl_EGlpNumIsLess (d_max, pinf->pmpinfo.infeas[j]))
      {
        dbl_EGlpNumCopy (d_max, pinf->pmpinfo.infeas[j]);
        pr->eindex = pinf->pmpinfo.bucket[j];
      }
    }
  }

  if (pr->eindex < 0)
    pr->price_stat = PRICE_OPTIMAL;
  else
  {
    if (phase == PRIMAL_PHASEI)
      dbl_EGlpNumCopy (d_e, lp->pIdz[pr->eindex]);
    else
      dbl_EGlpNumCopy (d_e, lp->dz[pr->eindex]);
    vs = lp->vstat[lp->nbaz[pr->eindex]];

    pr->price_stat = PRICE_NONOPTIMAL;
    if (vs == STAT_UPPER || (vs == STAT_ZERO && dbl_EGlpNumIsLess (*ftol, d_e)))
      pr->dir = VDECREASE;
    else
      pr->dir = VINCREASE;
  }
  dbl_EGlpNumClearVar (d_e);
  dbl_EGlpNumClearVar (d_max);
}

static void dbl_update_p_scaleinf (
  dbl_price_info * const p,
  dbl_heap * const h,
  int const i,
  double inf,
  int const prule)
{
  if (!dbl_EGlpNumIsNeqqZero (inf))
  {
    dbl_EGlpNumZero (p->p_scaleinf[i]);
    if (h->hexist != 0 && h->loc[i] != -1)
      dbl_ILLheap_delete (h, i);
  }
  else
  {
    if (prule == QS_PRICE_DDANTZIG)
      dbl_EGlpNumCopy (p->p_scaleinf[i], inf);
    else if (prule == QS_PRICE_DSTEEP)
      dbl_EGlpNumCopySqrOver (p->p_scaleinf[i], inf, p->dsinfo.norms[i]);
    else if (prule == QS_PRICE_DDEVEX)
      dbl_EGlpNumCopySqrOver (p->p_scaleinf[i], inf, p->ddinfo.norms[i]);

    if (h->hexist != 0)
    {
      if (h->loc[i] == -1)
        dbl_ILLheap_insert (h, i);
      else
        dbl_ILLheap_modify (h, i);
    }
  }
}

static void dbl_compute_primalI_inf (
  dbl_lpinfo * const lp,
  int const i,
  double * const inf)
{
  int const col = lp->baz[i];
  double*x = &(lp->xbz[i]);
  double*l = &(lp->lz[col]);
  double*u = &(lp->uz[col]);
  dbl_EGlpNumZero (*inf);
  double*ftol = &(lp->tol->ip_tol);

  if (dbl_EGlpNumIsLess (*ftol, *x) && dbl_EGlpNumIsNeqq (*u, dbl_INFTY))
    dbl_EGlpNumCopy (*inf, *x);
  else if (dbl_EGlpNumIsNeqq (*l, dbl_NINFTY) && dbl_EGlpNumIsSumLess (*x, *ftol,dbl_zeroLpNum))
    dbl_EGlpNumCopy (*inf, *x);
}

static void dbl_compute_primalII_inf (
  dbl_lpinfo * const lp,
  int const i,
  double * const inf)
{
  int const col = lp->baz[i];
  double*x = &(lp->xbz[i]);
  double*l = &(lp->lz[col]);
  double*u = &(lp->uz[col]);
  dbl_EGlpNumZero (*inf);
  double*ftol = &(lp->tol->pfeas_tol);

  if (dbl_EGlpNumIsNeqq (*u, dbl_INFTY) && dbl_EGlpNumIsSumLess (*u, *ftol, *x))
    dbl_EGlpNumCopyDiff (*inf, *x, *u);
  else if (dbl_EGlpNumIsNeqq (*l, dbl_NINFTY) && dbl_EGlpNumIsSumLess (*x, *ftol, *l))
    dbl_EGlpNumCopyDiff (*inf, *l, *x);
}

void dbl_ILLprice_compute_primal_inf (
  dbl_lpinfo * const lp,
  dbl_price_info * const p,
  int *const ix,
  int const icnt,
  int const phase)
{
  int i;
  int price;
  double inf;
  dbl_heap *h = &(p->h);

  price = (phase == DUAL_PHASEI) ? p->dI_price : p->dII_price;
  dbl_EGlpNumInitVar (inf);
  dbl_EGlpNumZero (inf);

  if (phase == DUAL_PHASEI)
  {
    if (ix == NULL)
      for (i = 0; i < lp->nrows; i++)
      {
        dbl_compute_primalI_inf (lp, i, &inf);
        dbl_update_p_scaleinf (p, h, i, inf, price);
      }
    else
      for (i = 0; i < icnt; i++)
      {
        dbl_compute_primalI_inf (lp, ix[i], &inf);
        dbl_update_p_scaleinf (p, h, ix[i], inf, price);
      }
  }
  else if (phase == DUAL_PHASEII)
  {
    if (ix == NULL)
      for (i = 0; i < lp->nrows; i++)
      {
        dbl_compute_primalII_inf (lp, i, &inf);
        dbl_update_p_scaleinf (p, h, i, inf, price);
      }
    else
      for (i = 0; i < icnt; i++)
      {
        dbl_compute_primalII_inf (lp, ix[i], &inf);
        dbl_update_p_scaleinf (p, h, ix[i], inf, price);
      }
  }
  dbl_EGlpNumClearVar (inf);
}

void dbl_ILLprice_dual (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  int const phase,
  dbl_price_res * const pr)
{
  int i;
  double p_max;
  double ubound;
  double*ftol = &(lp->tol->pfeas_tol);
  dbl_heap *const h = &(pinf->h);

  dbl_EGlpNumInitVar (p_max);
  dbl_EGlpNumInitVar (ubound);
  pr->lindex = -1;
  dbl_EGlpNumZero(p_max);

#if USEHEAP > 0
  dbl_ILLprice_test_for_heap (lp, pinf, lp->nrows, pinf->p_scaleinf, DUAL_SIMPLEX,
                          1);
#endif

  if (pinf->d_strategy == COMPLETE_PRICING)
  {
    if (h->hexist)
    {
      pr->lindex = dbl_ILLheap_findmin (h);
      if (pr->lindex != -1)
        dbl_ILLheap_delete (h, pr->lindex);
    }
    else
    {
      for (i = 0; i < lp->nrows; i++)
      {
        if (dbl_EGlpNumIsLess (p_max, pinf->p_scaleinf[i]))
        {
          dbl_EGlpNumCopy (p_max, pinf->p_scaleinf[i]);
          pr->lindex = i;
        }
      }
    }
  }
  else if (pinf->d_strategy == MULTI_PART_PRICING)
  {
    for (i = 0; i < pinf->dmpinfo.bsize; i++)
    {
      if (dbl_EGlpNumIsLess (p_max, pinf->dmpinfo.infeas[i]))
      {
        dbl_EGlpNumCopy (p_max, pinf->dmpinfo.infeas[i]);
        pr->lindex = pinf->dmpinfo.bucket[i];
      }
    }
  }

  if (pr->lindex < 0)
    pr->price_stat = PRICE_OPTIMAL;
  else
  {
    pr->price_stat = NONOPTIMAL;

    if (dbl_EGlpNumIsNeqq (lp->uz[lp->baz[pr->lindex]], dbl_INFTY))
    {
      if (phase == DUAL_PHASEI)
        dbl_EGlpNumZero(ubound);
      else
        dbl_EGlpNumCopy(ubound,lp->uz[lp->baz[pr->lindex]]);
      if (dbl_EGlpNumIsSumLess (*ftol, ubound, lp->xbz[pr->lindex]))
        pr->lvstat = STAT_UPPER;
      else
        pr->lvstat = STAT_LOWER;
    }
    else
      pr->lvstat = STAT_LOWER;
  }
  dbl_EGlpNumClearVar (p_max);
  dbl_EGlpNumClearVar (ubound);
}

int dbl_ILLprice_get_rownorms (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  double * const rnorms)
{
  int rval = 0;
  int i;

  if (pinf->dsinfo.norms == NULL)
  {
    rval = dbl_ILLprice_build_dsteep_norms (lp, &(pinf->dsinfo));
    CHECKRVALG(rval,CLEANUP);
  }
  for (i = 0; i < lp->nrows; i++)
    dbl_EGlpNumCopy (rnorms[i], pinf->dsinfo.norms[i]);

CLEANUP:
  if (rval)
    dbl_EGlpNumFreeArray (pinf->dsinfo.norms);

  return rval;
}

int dbl_ILLprice_get_colnorms (
  dbl_lpinfo * const lp,
  dbl_price_info * const pinf,
  double * const cnorms)
{
  int rval = 0;
  int i, j;

  if (pinf->psinfo.norms == NULL)
  {
    rval = dbl_ILLprice_build_psteep_norms (lp, &(pinf->psinfo));
    CHECKRVALG(rval,CLEANUP);
  }
  for (i = 0; i < lp->nrows; i++)
    dbl_EGlpNumZero (cnorms[lp->baz[i]]);
  for (j = 0; j < lp->nnbasic; j++)
    dbl_EGlpNumCopy (cnorms[lp->nbaz[j]], pinf->psinfo.norms[j]);

CLEANUP:
  if (rval)
    dbl_EGlpNumFreeArray (pinf->psinfo.norms);

  return rval;
}

int dbl_ILLprice_get_newnorms (
  dbl_lpinfo * const lp,
  int const nelems,
  double * const norms,
  int *const matcnt,
  int *const matbeg,
  int *const matind,
  double * const matval,
  int const option)
{
  int i, j;
  int rval = 0;
  dbl_svector a;
  dbl_svector y;

  dbl_ILLsvector_init (&y);
  rval = dbl_ILLsvector_alloc (&y, lp->nrows);
  CHECKRVALG(rval,CLEANUP);

  for (j = 0; j < nelems; j++)
  {
    a.nzcnt = matcnt[j];
    a.indx = &(matind[matbeg[j]]);
    a.coef = &(matval[matbeg[j]]);

    if (option == COLUMN_SOLVE)
      dbl_ILLbasis_column_solve (lp, &a, &y);
    else
      dbl_ILLbasis_row_solve (lp, &a, &y);

    dbl_EGlpNumOne (norms[j]);
    for (i = 0; i < y.nzcnt; i++)
      dbl_EGlpNumAddInnProdTo (norms[j], y.coef[i], y.coef[i]);
  }

CLEANUP:
  dbl_ILLsvector_free (&y);
  EG_RETURN(rval);
}

int dbl_ILLprice_get_new_rownorms (
  dbl_lpinfo * const lp,
  int const newrows,
  double * const rnorms,
  int *const rmatcnt,
  int *const rmatbeg,
  int *const rmatind,
  double * const rmatval)
{
  return dbl_ILLprice_get_newnorms (lp, newrows, rnorms, rmatcnt, rmatbeg, rmatind,
                                rmatval, ROW_SOLVE);
}

int dbl_ILLprice_get_new_colnorms (
  dbl_lpinfo * const lp,
  int const newrows,
  double * const rnorms,
  int *const matcnt,
  int *const matbeg,
  int *const matind,
  double * const matval)
{
  return dbl_ILLprice_get_newnorms (lp, newrows, rnorms, matcnt, matbeg, matind,
                                matval, COLUMN_SOLVE);
}

int dbl_ILLprice_load_rownorms (
  dbl_lpinfo * const lp,
  double * const rnorms,
  dbl_price_info * const pinf)
{
  int i;
  int rval = 0;

  dbl_EGlpNumFreeArray (pinf->dsinfo.norms);
  pinf->dsinfo.norms = dbl_EGlpNumAllocArray (lp->nrows);

  for (i = 0; i < lp->nrows; i++)
  {
    dbl_EGlpNumCopy (pinf->dsinfo.norms[i], rnorms[i]);
    if (dbl_EGlpNumIsLess (pinf->dsinfo.norms[i], dbl_PARAM_MIN_DNORM))
      dbl_EGlpNumCopy (pinf->dsinfo.norms[i], dbl_PARAM_MIN_DNORM);
  }

  EG_RETURN(rval);
}

int dbl_ILLprice_load_colnorms (
  dbl_lpinfo * const lp,
  double * const cnorms,
  dbl_price_info * const pinf)
{
  int j;
  int rval = 0;

  dbl_EGlpNumFreeArray (pinf->psinfo.norms);
  pinf->psinfo.norms = dbl_EGlpNumAllocArray (lp->nnbasic);

  for (j = 0; j < lp->nnbasic; j++)
  {
    dbl_EGlpNumCopy (pinf->psinfo.norms[j], cnorms[lp->nbaz[j]]);
    if (dbl_EGlpNumIsLess (pinf->psinfo.norms[j], dbl_oneLpNum))
      dbl_EGlpNumOne (pinf->psinfo.norms[j]);
  }

  EG_RETURN(rval);
}

#if dbl_PRICE_DEBUG > 0
void dbl_test_dsteep_norms (
  dbl_lpinfo * lp,
  dbl_price_info * p)
{
  int i, errn = 0;
  double *pn = dbl_EGlpNumAllocArray(lp->nrows);
  double err, diff;
  dbl_EGlpNumZero (err);

  dbl_EGlpNumInitVar (err);
  dbl_EGlpNumInitVar (diff);

  dbl_ILLprice_get_dsteep_norms (lp, lp->yjz.nzcnt, lp->yjz.indx, pn);
  for (i = 0; i < lp->yjz.nzcnt; i++)
  {
    dbl_EGlpNumCopyDiff (diff, pn[i], p->dsinfo.norms[lp->yjz.indx[i]]);
    dbl_EGlpNumCopyAbs(diff,diff);
    if (dbl_EGlpNumIsLess (dbl_PFEAS_TOLER, diff))
    {
      errn++;
      dbl_EGlpNumAddTo (err, diff);
      dbl_EGlpNumCopy (p->dsinfo.norms[lp->yjz.indx[i]], pn[i]);
    }
  }
  if (errn)
    printf ("%d: dnorm errn = %d, err = %.6f\n", lp->cnts->tot_iter, errn,
            dbl_EGlpNumToLf (err));
  dbl_EGlpNumFreeArray (pn);
  dbl_EGlpNumClearVar (diff);
  dbl_EGlpNumClearVar (err);
}
#endif

---