CoinOslC.h

Go to the documentation of this file.

/* $Id: CoinOslC.h 2259 2020-01-16 13:43:43Z stefan $ */
#ifndef COIN_OSL_C_INCLUDE
/*
  Copyright (C) 1987, 2009, International Business Machines Corporation
  and others.  All Rights Reserved.

  This code is licensed under the terms of the Eclipse Public License (EPL).
*/
#define COIN_OSL_C_INCLUDE

#ifndef CLP_OSL
#define CLP_OSL 0
#endif
#define C_EKK_GO_SPARSE 200

#ifdef HAVE_ENDIAN_H
#include <endian.h>
#if __BYTE_ORDER == __LITTLE_ENDIAN
#define INTEL
#endif
#endif

#include <math.h>
#include <string.h>
#include <stdio.h>
#include <stdlib.h>

#define SPARSE_UPDATE
#define NO_SHIFT
#include "CoinHelperFunctions.hpp"

#include <stddef.h>
#ifdef __cplusplus
extern "C" {
#endif

int c_ekkbtrn(register const EKKfactinfo *fact,
  double *dwork1,
  int *mpt, int first_nonzero);
int c_ekkbtrn_ipivrw(register const EKKfactinfo *fact,
  double *dwork1,
  int *mpt, int ipivrw, int *spare);

int c_ekketsj(register /*const*/ EKKfactinfo *fact,
  double *dwork1,
  int *mpt2, double dalpha, int orig_nincol,
  int npivot, int *nuspikp,
  const int ipivrw, int *spare);
int c_ekkftrn(register const EKKfactinfo *fact,
  double *dwork1,
  double *dpermu, int *mpt, int numberNonZero);

int c_ekkftrn_ft(register EKKfactinfo *fact,
  double *dwork1, int *mpt, int *nincolp);
void c_ekkftrn2(register EKKfactinfo *fact, double *dwork1,
  double *dpermu1, int *mpt1, int *nincolp,
  double *dwork1_ft, int *mpt_ft, int *nincolp_ft);

int c_ekklfct(register EKKfactinfo *fact);
int c_ekkslcf(register const EKKfactinfo *fact);
inline void c_ekkscpy(int n, const int *marr1, int *marr2)
{
  CoinMemcpyN(marr1, n, marr2);
}
inline void c_ekkdcpy(int n, const double *marr1, double *marr2)
{
  CoinMemcpyN(marr1, n, marr2);
}
int c_ekk_IsSet(const int *array, int bit);
void c_ekk_Set(int *array, int bit);
void c_ekk_Unset(int *array, int bit);

void c_ekkzero(int length, int n, void *array);
inline void c_ekkdzero(int n, double *marray)
{
  CoinZeroN(marray, n);
}
inline void c_ekkizero(int n, int *marray)
{
  CoinZeroN(marray, n);
}
inline void c_ekkczero(int n, char *marray)
{
  CoinZeroN(marray, n);
}
#ifdef __cplusplus
}
#endif

#define c_ekkscpy_0_1(s, ival, array) CoinFillN(array, s, ival)
#define c_ekks1cpy(n, marr1, marr2) CoinMemcpyN(marr1, n, marr2)
void clp_setup_pointers(EKKfactinfo *fact);
void clp_memory(int type);
double *clp_double(int number_entries);
int *clp_int(int number_entries);
void *clp_malloc(int number_entries);
void clp_free(void *oldArray);

#define SLACK_VALUE -1.0
#define C_EKK_REMOVE_LINK(hpiv, hin, link, ipivot) \
  {                                                \
    int ipre = link[ipivot].pre;                   \
    int isuc = link[ipivot].suc;                   \
    if (ipre > 0) {                                \
      link[ipre].suc = isuc;                       \
    }                                              \
    if (ipre <= 0) {                               \
      hpiv[hin[ipivot]] = isuc;                    \
    }                                              \
    if (isuc > 0) {                                \
      link[isuc].pre = ipre;                       \
    }                                              \
  }

#define C_EKK_ADD_LINK(hpiv, nzi, link, npr) \
  {                                          \
    int ifiri = hpiv[nzi];                   \
    hpiv[nzi] = npr;                         \
    link[npr].suc = ifiri;                   \
    link[npr].pre = 0;                       \
    if (ifiri != 0) {                        \
      link[ifiri].pre = npr;                 \
    }                                        \
  }
#include <assert.h>
#ifdef NO_SHIFT

#define SHIFT_INDEX(limit) (limit)
#define UNSHIFT_INDEX(limit) (limit)
#define SHIFT_REF(arr, ind) (arr)[ind]

#else

#define SHIFT_INDEX(limit) ((limit) << 3)
#define UNSHIFT_INDEX(limit) ((unsigned int)(limit) >> 3)
#define SHIFT_REF(arr, ind) (*(double *)((char *)(arr) + (ind)))

#endif

#ifdef INTEL
#define NOT_ZERO(x) (((*((reinterpret_cast< unsigned char * >(&x)) + 7)) & 0x7F) != 0)
#else
#define NOT_ZERO(x) ((x) != 0.0)
#endif

#define SWAP(type, _x, _y) \
  {                        \
    type _tmp = (_x);      \
    (_x) = (_y);           \
    (_y) = _tmp;           \
  }

#define UNROLL_LOOP_BODY1(code) \
  {                             \
    {                           \
      code                      \
    }                           \
  }
#define UNROLL_LOOP_BODY2(code) \
  {                             \
    { code } {                  \
      code                      \
    }                           \
  }
#define UNROLL_LOOP_BODY4(code)  \
  {                              \
    { code } { code } { code } { \
      code                       \
    }                            \
  }
#endif
#ifdef COIN_OSL_CMFC
/*     Return codes in IRTCOD/IRTCOD are */
/*     4: numerical problems */
/*     5: not enough space in row file */
/*     6: not enough space in column file */
/*    23: system error at label 320 */
{
#if 1
  int *hcoli = fact->xecadr;
  double *dluval = fact->xeeadr;
  double *dvalpv = fact->kw3adr;
  int *mrstrt = fact->xrsadr;
  int *hrowi = fact->xeradr;
  int *mcstrt = fact->xcsadr;
  int *hinrow = fact->xrnadr;
  int *hincol = fact->xcnadr;
  int *hpivro = fact->krpadr;
  int *hpivco = fact->kcpadr;
#endif
  int nnentl = fact->nnentl;
  int nnentu = fact->nnentu;
  int kmxeta = fact->kmxeta;
  int xnewro = *xnewrop;
  int ncompactions = *ncompactionsp;

  MACTION_T *maction = reinterpret_cast< MACTION_T * >(maction_void);

  int i, j, k;
  double d1;
  int j1, j2;
  int jj, kk, kr, nz, jj1, jj2, kce, kcs, kqq, npr;
  int fill, naft;
  int enpr;
  int nres, npre;
  int knpr, irow, iadd32, ibase;
  double pivot;
  int count, nznpr;
  int nlast, epivr1;
  int kipis;
  double dpivx;
  int kipie, kcpiv, knprs, knpre;
  bool cancel;
  double multip, elemnt;
  int ipivot, jpivot, epivro, epivco, lstart, nfirst;
  int nzpivj, kfill, kstart;
  int nmove, ileft;
#ifndef C_EKKCMFY
  int iput, nspare;
  int noRoomForDense = 0;
  int if_sparse_update = fact->if_sparse_update;
  int ifdens = 0;
#endif
  int irtcod = 0;
  const int nrow = fact->nrow;

  /* Parameter adjustments */
  --maction;

  /* Function Body */
  lstart = nnetas - nnentl + 1;
  for (i = lstart; i <= nnetas; ++i) {
    hrowi[i] = SHIFT_INDEX(hcoli[i]);
  }

  for (i = 1; i <= nrow; ++i) {
    maction[i] = 0;
    mwork[i].pre = i - 1;
    mwork[i].suc = i + 1;
  }

  iadd32 = 0;
  nlast = nrow;
  nfirst = 1;
  mwork[1].pre = nrow;
  mwork[nrow].suc = 1;

  for (count = 1; count <= nrow; ++count) {

    /* Pick column singletons */
    if (!(hpivco[1] <= 0)) {
      int small_pivot = c_ekkcsin(fact,
        rlink, clink,
        nsingp);

      if (small_pivot) {
        irtcod = 7; /* pivot too small */
        if (fact->invok >= 0) {
          goto L1050;
        }
      }
      if (fact->npivots >= nrow) {
        goto L1050;
      }
    }

    /* Pick row singletons */
    if (!(hpivro[1] <= 0)) {
      irtcod = c_ekkrsin(fact,
        rlink, clink,
        mwork, nfirst,
        nsingp,

        &xnewco, &xnewro,
        &nnentu,
        &kmxeta, &ncompactions,
        &nnentl);
      if (irtcod != 0) {
        if (irtcod < 0 || fact->invok >= 0) {
          /* -5 */
          goto L1050;
        }
        /* ASSERT:  irtcod == 7 - pivot too small */
        /* why don't we return with an error? */
      }
      if (fact->npivots >= nrow) {
        goto L1050;
      }
      lstart = nnetas - nnentl + 1;
    }

    /* Find a pivot element */
    irtcod = c_ekkfpvt(fact,
      rlink, clink,
      nsingp, xrejctp, &ipivot, &jpivot);
    if (irtcod != 0) {
      /* irtcod == 10 */
      goto L1050;
    }
    /*        Update list structures and prepare for numerical phase */
    c_ekkprpv(fact, rlink, clink,
      *xrejctp, ipivot, jpivot);

    epivco = hincol[jpivot];
    ++fact->xnetal;
    mcstrt[fact->xnetal] = lstart - 1;
    hpivco[fact->xnetal] = ipivot;
    epivro = hinrow[ipivot];
    epivr1 = epivro - 1;
    kipis = mrstrt[ipivot];
    pivot = dluval[kipis];
    dpivx = 1. / pivot;
    kipie = kipis + epivr1;
    ++kipis;
#ifndef C_EKKCMFY
    {
      double size = nrow - fact->npivots;
      if (size > GO_DENSE && (nnentu - fact->nuspike) * GO_DENSE_RATIO > size * size) {
        /* say going to dense coding */
        if (*nsingp == 0) {
          ifdens = 1;
        }
      }
    }
#endif
    /* copy the pivot row entries into dvalpv */
    /* the maction array tells us the index into dvalpv for a given row */
    /* the alternative would be using a large array of doubles */
    for (k = kipis; k <= kipie; ++k) {
      irow = hcoli[k];
      dvalpv[k - kipis + 1] = dluval[k];
      maction[irow] = static_cast< MACTION_T >(k - kipis + 1);
    }

    /* Loop over nonzeros in pivot column */
    kcpiv = mcstrt[jpivot] - 1;
    for (nzpivj = 1; nzpivj <= epivco; ++nzpivj) {
      ++kcpiv;
      npr = hrowi[kcpiv];
      hrowi[kcpiv] = 0; /* zero out for possible compaction later on */

      --hincol[jpivot];

      ++mcstrt[jpivot];
      /* loop invariant:  kcpiv == mcstrt[jpivot] - 1 */

      --hinrow[npr];
      enpr = hinrow[npr];
      knprs = mrstrt[npr];
      knpre = knprs + enpr;

      /* Search for element to be eliminated */
      knpr = knprs;
      while (1) {
        UNROLL_LOOP_BODY4({
          if (jpivot == hcoli[knpr]) {
            break;
          }
          knpr++;
        });
      }

      multip = -dluval[knpr] * dpivx;

      /* swap last entry with pivot */
      dluval[knpr] = dluval[knpre];
      hcoli[knpr] = hcoli[knpre];
      --knpre;

#if 1
      /* MONSTER_UNROLLED_CODE - see below */
      kfill = epivr1 - (knpre - knprs + 1);
      nres = ((knpre - knprs + 1) & 1) + knprs;
      cancel = false;
      d1 = 1e33;
      j1 = hcoli[nres];

      if (nres != knprs) {
        j = hcoli[knprs];
        if (maction[j] == 0) {
          ++kfill;
        } else {
          jj = maction[j];
          maction[j] = static_cast< MACTION_T >(-maction[j]);
          dluval[knprs] += multip * dvalpv[jj];
          d1 = fabs(dluval[knprs]);
        }
      }
      j2 = hcoli[nres + 1];
      jj1 = maction[j1];
      for (kr = nres; kr < knpre; kr += 2) {
        jj2 = maction[j2];
        if (jj1 == 0) {
          ++kfill;
        } else {
          maction[j1] = static_cast< MACTION_T >(-maction[j1]);
          dluval[kr] += multip * dvalpv[jj1];
          cancel = cancel || !(fact->zeroTolerance < d1);
          d1 = fabs(dluval[kr]);
        }
        j1 = hcoli[kr + 2];
        if (jj2 == 0) {
          ++kfill;
        } else {
          maction[j2] = static_cast< MACTION_T >(-maction[j2]);
          dluval[kr + 1] += multip * dvalpv[jj2];
          cancel = cancel || !(fact->zeroTolerance < d1);
          d1 = fabs(dluval[kr + 1]);
        }
        jj1 = maction[j1];
        j2 = hcoli[kr + 3];
      }
      cancel = cancel || !(fact->zeroTolerance < d1);
#else
      /*
       * This is apparently what the above code does.
       * In addition to being unrolled, the assignments to j[12] and jj[12]
       * are shifted so that the result of dereferencing maction doesn't
       * have to be used immediately afterwards for the branch test.
       * This would would cause a pipeline delay.  (The apparent dereference
       * of hcoli will be removed by the compiler using strength reduction).
       *
       * loop through the entries in the row being processed,
       * flipping the sign of the maction entries as we go along.
       * Afterwards, we look for positive entries to see what pivot
       * row entries will cause fill-in.  We count the number of fill-ins, too.
       * "cancel" says if the result of combining the pivot row with this one
       * causes an entry to get too small; if so, we discard those entries.
       */
      kfill = epivr1 - (knpre - knprs + 1);
      cancel = false;

      for (kr = knprs; kr <= knpre; kr++) {
        j1 = hcoli[kr];
        jj1 = maction[j1];
        if ((jj1 == 0)) {
          /* no entry - this pivot column entry will have to be added */
          ++kfill;
        } else {
          /* there is an entry for this column in the pivot row */
          maction[j1] = -maction[j1];
          dluval[kr] += multip * dvalpv[jj1];
          d1 = fabs(dluval[kr]);
          cancel = cancel || !(fact->zeroTolerance < d1);
        }
      }
#endif
      kstart = knpre;
      fill = kfill;

      if (cancel) {
        /* KSTART is used as a stack pointer for nonzeros in factored row */
        kstart = knprs - 1;
        for (kr = knprs; kr <= knpre; ++kr) {
          j = hcoli[kr];
          if (fabs(dluval[kr]) > fact->zeroTolerance) {
            ++kstart;
            dluval[kstart] = dluval[kr];
            hcoli[kstart] = j;
          } else {
            /* Remove element from column file */
            --nnentu;
            --hincol[j];
            --enpr;
            kcs = mcstrt[j];
            kce = kcs + hincol[j];
            for (kk = kcs; kk <= kce; ++kk) {
              if (hrowi[kk] == npr) {
                hrowi[kk] = hrowi[kce];
                hrowi[kce] = 0;
                break;
              }
            }
            /* ASSERT !(kk>kce) */
          }
        }
        knpre = kstart;
      }
      /* Fill contains an upper bound on the amount of fill-in */
      if (fill == 0) {
        for (k = kipis; k <= kipie; ++k) {
          maction[hcoli[k]] = static_cast< MACTION_T >(-maction[hcoli[k]]);
        }
      } else {
        naft = mwork[npr].suc;
        kqq = mrstrt[naft] - knpre - 1;

        if (fill > kqq) {
          /* Fill-in exceeds space left. Check if there is enough */
          /* space in row file for the new row. */
          nznpr = enpr + fill;
          if (!(xnewro + nznpr + 1 < lstart)) {
            if (!(nnentu + nznpr + 1 < lstart)) {
              irtcod = -5;
              goto L1050;
            }
            /* idea 1 is to compress every time xnewro increases by x thousand */
            /* idea 2 is to copy nucleus rows with a reasonable gap */
            /* then copy each row down when used */
            /* compressions would just be 1 remainder which eventually will */
            /* fit in cache */
            {
              int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
              kmxeta += xnewro - iput;
              xnewro = iput - 1;
              ++ncompactions;
            }

            kipis = mrstrt[ipivot] + 1;
            kipie = kipis + epivr1 - 1;
            knprs = mrstrt[npr];
          }

          /* I think this assignment should be inside the previous if-stmt */
          /* otherwise, it does nothing */
          /*assert(knpre == knprs + enpr - 1);*/
          knpre = knprs + enpr - 1;

          /*
           * copy this row to the end of the row file and adjust its links.
           * The links keep track of the order of rows in memory.
           * Rows are only moved from the middle all the way to the end.
           */
          if (npr != nlast) {
            npre = mwork[npr].pre;
            if (npr == nfirst) {
              nfirst = naft;
            }
            /*             take out of chain */
            mwork[naft].pre = npre;
            mwork[npre].suc = naft;
            /*             and put in at end */
            mwork[nfirst].pre = npr;
            mwork[nlast].suc = npr;
            mwork[npr].pre = nlast;
            mwork[npr].suc = nfirst;
            nlast = npr;
            kstart = xnewro;
            mrstrt[npr] = kstart + 1;
            nmove = knpre - knprs + 1;
            ibase = kstart + 1 - knprs;
            for (kr = knprs; kr <= knpre; ++kr) {
              dluval[ibase + kr] = dluval[kr];
              hcoli[ibase + kr] = hcoli[kr];
            }
            kstart += nmove;
          } else {
            kstart = knpre;
          }

          /* extra space ? */
          /*
           * The mystery of iadd32.
           * This code assigns to xnewro, possibly using iadd32.
           * However, in that case xnewro is assigned to just after
           * the for-loop below, and there is no intervening reference.
           * Therefore, I believe that this code can be entirely eliminated;
           * it is the leftover of an interrupted or dropped experiment.
           * Presumably, this was trying to implement the ideas about
           * padding expressed above.
           */
          if (iadd32 != 0) {
            xnewro += iadd32;
          } else {
            if (kstart + (nrow << 1) + 100 < lstart) {
              ileft = ((nrow - fact->npivots + 32) & -32);
              if (kstart + ileft * ileft + 32 < lstart) {
                iadd32 = ileft;
                xnewro = CoinMax(kstart, xnewro);
                xnewro = (xnewro & -32) + ileft;
              } else {
                xnewro = ((kstart + 31) & -32);
              }
            } else {
              xnewro = kstart;
            }
          }

          hinrow[npr] = enpr;
        } else if (!(nnentu + kqq + 2 < lstart)) {
          irtcod = -5;
          goto L1050;
        }
        /* Scan pivot row again to generate fill in. */
        for (kr = kipis; kr <= kipie; ++kr) {
          j = hcoli[kr];
          jj = maction[j];
          if (jj > 0) {
            elemnt = multip * dvalpv[jj];
            if (fabs(elemnt) > fact->zeroTolerance) {
              ++kstart;
              dluval[kstart] = elemnt;
              //printf("pivot %d at %d col %d el %g\n",
              // npr,kstart,j,elemnt);
              hcoli[kstart] = j;
              ++nnentu;
              nz = hincol[j];
              kcs = mcstrt[j];
              kce = kcs + nz - 1;
              if (kce == xnewco) {
                if (xnewco + 1 >= lstart) {
                  if (xnewco + nz + 1 >= lstart) {
                    /*                  Compress column file */
                    if (nnentu + nz + 1 < lstart) {
                      xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
                      ++ncompactions;

                      kcpiv = mcstrt[jpivot] - 1;
                      kcs = mcstrt[j];
                      /*                  HINCOL MAY HAVE CHANGED? (JJHF) ??? */
                      nz = hincol[j];
                      kce = kcs + nz - 1;
                    } else {
                      irtcod = -5;
                      goto L1050;
                    }
                  }
                  /*              Copy column */
                  mcstrt[j] = xnewco + 1;
                  ibase = mcstrt[j] - kcs;
                  for (kk = kcs; kk <= kce; ++kk) {
                    hrowi[ibase + kk] = hrowi[kk];
                    hrowi[kk] = 0;
                  }
                  kce = xnewco + kce - kcs + 1;
                  xnewco = kce + 1;
                } else {
                  ++xnewco;
                }
              } else if (hrowi[kce + 1] != 0) {
                /* here we use the fact that hrowi entries not "in use" are zeroed */
                if (xnewco + nz + 1 >= lstart) {
                  /* Compress column file */
                  if (nnentu + nz + 1 < lstart) {
                    xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
                    ++ncompactions;

                    kcpiv = mcstrt[jpivot] - 1;
                    kcs = mcstrt[j];
                    /*                  HINCOL MAY HAVE CHANGED? (JJHF) ??? */
                    nz = hincol[j];
                    kce = kcs + nz - 1;
                  } else {
                    irtcod = -5;
                    goto L1050;
                  }
                }
                /* move the column to the end of the column file */
                mcstrt[j] = xnewco + 1;
                ibase = mcstrt[j] - kcs;
                for (kk = kcs; kk <= kce; ++kk) {
                  hrowi[ibase + kk] = hrowi[kk];
                  hrowi[kk] = 0;
                }
                kce = xnewco + kce - kcs + 1;
                xnewco = kce + 1;
              }
              /* store element */
              hrowi[kce + 1] = npr;
              hincol[j] = nz + 1;
            }
          } else {
            maction[j] = static_cast< MACTION_T >(-maction[j]);
          }
        }
        if (fill > kqq) {
          xnewro = kstart;
        }
      }
      hinrow[npr] = kstart - mrstrt[npr] + 1;
      /* Check if row or column file needs compression */
      if (!(xnewco + 1 < lstart)) {
        xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
        ++ncompactions;

        kcpiv = mcstrt[jpivot] - 1;
      }
      if (!(xnewro + 1 < lstart)) {
        int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
        kmxeta += xnewro - iput;
        xnewro = iput - 1;
        ++ncompactions;

        kipis = mrstrt[ipivot] + 1;
        kipie = kipis + epivr1 - 1;
      }
      /* Store elementary row transformation */
      ++nnentl;
      --nnentu;
      --lstart;
      dluval[lstart] = multip;

      hrowi[lstart] = SHIFT_INDEX(npr);
#define INLINE_AFPV 3
      /* We could do this while computing values but
         it makes it much more complex.  At least we should get
         reasonable cache behavior by doing it each row */
#if INLINE_AFPV
      {
        int j;
        int nel, krs;
        int koff;
        int *index;
        double *els;
        nel = hinrow[npr];
        krs = mrstrt[npr];
        index = &hcoli[krs];
        els = &dluval[krs];
#if INLINE_AFPV < 3
#if INLINE_AFPV == 1
        double maxaij = 0.0;
        koff = 0;
        j = 0;
        while (j < nel) {
          double d = fabs(els[j]);
          if (maxaij < d) {
            maxaij = d;
            koff = j;
          }
          j++;
        }
#else
        assert(nel);
        koff = 0;
        double maxaij = fabs(els[0]);
        for (j = 1; j < nel; j++) {
          double d = fabs(els[j]);
          if (maxaij < d) {
            maxaij = d;
            koff = j;
          }
        }
#endif
#else
        double maxaij = 0.0;
        koff = 0;
        j = 0;
        if ((nel & 1) != 0) {
          maxaij = fabs(els[0]);
          j = 1;
        }

        while (j < nel) {
          UNROLL_LOOP_BODY2({
            double d = fabs(els[j]);
            if (maxaij < d) {
              maxaij = d;
              koff = j;
            }
            j++;
          });
        }
#endif
        SWAP(int, index[koff], index[0]);
        SWAP(double, els[koff], els[0]);
      }
#endif

      {
        int nzi = hinrow[npr];
        if (nzi > 0) {
          C_EKK_ADD_LINK(hpivro, nzi, rlink, npr);
        }
      }
    }

    /* after pivot move biggest to first in each row */
#if INLINE_AFPV == 0
    int nn = mcstrt[fact->xnetal] - lstart + 1;
    c_ekkafpv(hrowi + lstart, hcoli, dluval, mrstrt, hinrow, nn);
#endif

    /* Restore work array */
    for (k = kipis; k <= kipie; ++k) {
      maction[hcoli[k]] = 0;
    }

    if (*xrejctp > 0) {
      for (k = kipis; k <= kipie; ++k) {
        int j = hcoli[k];
        int nzj = hincol[j];
        if (!(nzj <= 0) && !((clink[j].pre > nrow && nzj != 1))) {
          C_EKK_ADD_LINK(hpivco, nzj, clink, j);
        }
      }
    } else {
      for (k = kipis; k <= kipie; ++k) {
        int j = hcoli[k];
        int nzj = hincol[j];
        if (!(nzj <= 0)) {
          C_EKK_ADD_LINK(hpivco, nzj, clink, j);
        }
      }
    }
    fact->nuspike += hinrow[ipivot];

    /* Go to dense coding if appropriate */
#ifndef C_EKKCMFY
    if (ifdens != 0) {
      int ndense = nrow - fact->npivots;
      if (!(xnewro + ndense * ndense >= lstart)) {

        /* set up sort order in MACTION */
        c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
        iput = 0;
        for (i = 1; i <= nrow; ++i) {
          if (clink[i].pre >= 0) {
            ++iput;
            maction[i] = static_cast< short int >(iput);
          }
        }
        /* and get number spare needed */
        nspare = 0;
        for (i = 1; i <= nrow; ++i) {
          if (rlink[i].pre >= 0) {
            nspare = nspare + ndense - hinrow[i];
          }
        }
        if (iput != nrow - fact->npivots) {
          /* must be singular */
          c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
        } else {
          /* pack down then back up */
          int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
          kmxeta += xnewro - iput;
          xnewro = iput - 1;
          ++ncompactions;

          --ncompactions;
          if (xnewro + nspare + ndense * ndense >= lstart) {
            c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
          } else {
            xnewro += nspare;
            c_ekkrwct(fact, dluval, hcoli, mrstrt, hinrow, mwork,
              rlink, maction, dvalpv,
              nlast, xnewro);
            kmxeta += xnewro;
            if (nnentu + nnentl > nrow * 5 && (ndense * ndense) > (nnentu + nnentl) >> 2 && !if_sparse_update) {
              fact->ndenuc = ndense;
            }
            irtcod = c_ekkcmfd(fact,
              (reinterpret_cast< int * >(dvalpv) + 1),
              rlink, clink,
              (reinterpret_cast< int * >(maction + 1)) + 1,
              nnetas,
              &nnentl, &nnentu,
              nsingp);
            /* irtcod == 0 || irtcod == 10 */
            /* 10 == found 0.0 pivot */
            goto L1050;
          }
        }
      } else {
        /* say not enough room */
        /*printf("no room %d\n",ndense);*/
        if (1) {
          /* return and increase size of etas if possible */
          if (!noRoomForDense) {
            int etasize = CoinMax(4 * fact->nnentu + (nnetas - fact->nnentl) + 1000, fact->eta_size);
            noRoomForDense = ndense;
            fact->eta_size = CoinMin(static_cast< int >(1.2 * fact->eta_size), etasize);
            if (fact->maxNNetas > 0 && fact->eta_size > fact->maxNNetas) {
              fact->eta_size = fact->maxNNetas;
            }
          }
        }
      }
    }
#endif /* C_EKKCMFY */
  }

L1050 : {
  int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
  kmxeta += xnewro - iput;
  xnewro = iput - 1;
  ++ncompactions;
}

  nnentu = xnewro;
  /* save order of row copy for c_ekkshfv */
  mwork[nrow + 1].pre = nfirst;
  mwork[nrow + 1].suc = nlast;

  fact->nnentl = nnentl;
  fact->nnentu = nnentu;
  fact->kmxeta = kmxeta;
  *xnewrop = xnewro;
  *ncompactionsp = ncompactions;

  return (irtcod);
} /* c_ekkcmfc */
#endif

/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
*/