Ignore:
Timestamp:
Jan 6, 2019 2:43:06 PM (3 months ago)
Author:
unxusr
Message:

formatting

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/ClpNetworkBasis.cpp

    r2271 r2385  
    1111#include "CoinIndexedVector.hpp"
    1212
    13 
    1413//#############################################################################
    1514// Constructors / Destructor / Assignment
     
    1918// Default Constructor
    2019//-------------------------------------------------------------------
    21 ClpNetworkBasis::ClpNetworkBasis ()
     20ClpNetworkBasis::ClpNetworkBasis()
    2221{
    2322#ifndef COIN_FAST_CODE
    24      slackValue_ = -1.0;
     23  slackValue_ = -1.0;
    2524#endif
    26      numberRows_ = 0;
    27      numberColumns_ = 0;
    28      parent_ = NULL;
    29      descendant_ = NULL;
    30      pivot_ = NULL;
    31      rightSibling_ = NULL;
    32      leftSibling_ = NULL;
    33      sign_ = NULL;
    34      stack_ = NULL;
    35      permute_ = NULL;
    36      permuteBack_ = NULL;
    37      stack2_ = NULL;
    38      depth_ = NULL;
    39      mark_ = NULL;
    40      model_ = NULL;
     25  numberRows_ = 0;
     26  numberColumns_ = 0;
     27  parent_ = NULL;
     28  descendant_ = NULL;
     29  pivot_ = NULL;
     30  rightSibling_ = NULL;
     31  leftSibling_ = NULL;
     32  sign_ = NULL;
     33  stack_ = NULL;
     34  permute_ = NULL;
     35  permuteBack_ = NULL;
     36  stack2_ = NULL;
     37  depth_ = NULL;
     38  mark_ = NULL;
     39  model_ = NULL;
    4140}
    4241// Constructor from CoinFactorization
    43 ClpNetworkBasis::ClpNetworkBasis(const ClpSimplex * model,
    44                                  int numberRows, const CoinFactorizationDouble * pivotRegion,
    45                                  const int * permuteBack,
    46                                  const int * startColumn,
    47                                  const int * numberInColumn,
    48                                  const int * indexRow, const CoinFactorizationDouble * /*element*/)
     42ClpNetworkBasis::ClpNetworkBasis(const ClpSimplex *model,
     43  int numberRows, const CoinFactorizationDouble *pivotRegion,
     44  const int *permuteBack,
     45  const int *startColumn,
     46  const int *numberInColumn,
     47  const int *indexRow, const CoinFactorizationDouble * /*element*/)
    4948{
    5049#ifndef COIN_FAST_CODE
    51      slackValue_ = -1.0;
     50  slackValue_ = -1.0;
    5251#endif
    53      numberRows_ = numberRows;
    54      numberColumns_ = numberRows;
    55      parent_ = new int [ numberRows_+1];
    56      descendant_ = new int [ numberRows_+1];
    57      pivot_ = new int [ numberRows_+1];
    58      rightSibling_ = new int [ numberRows_+1];
    59      leftSibling_ = new int [ numberRows_+1];
    60      sign_ = new double [ numberRows_+1];
    61      stack_ = new int [ numberRows_+1];
    62      stack2_ = new int[numberRows_+1];
    63      depth_ = new int[numberRows_+1];
    64      mark_ = new char[numberRows_+1];
    65      permute_ = new int [numberRows_ + 1];
    66      permuteBack_ = new int [numberRows_ + 1];
    67      int i;
    68      for (i = 0; i < numberRows_ + 1; i++) {
    69           parent_[i] = -1;
    70           descendant_[i] = -1;
    71           pivot_[i] = -1;
    72           rightSibling_[i] = -1;
    73           leftSibling_[i] = -1;
    74           sign_[i] = -1.0;
    75           stack_[i] = -1;
    76           permute_[i] = i;
    77           permuteBack_[i] = i;
    78           stack2_[i] = -1;
    79           depth_[i] = -1;
    80           mark_[i] = 0;
    81      }
    82      mark_[numberRows_] = 1;
    83      // pivotColumnBack gives order of pivoting into basis
    84      // so pivotColumnback[0] is first slack in basis and
    85      // it pivots on row permuteBack[0]
    86      // a known root is given by permuteBack[numberRows_-1]
    87      for (i = 0; i < numberRows_; i++) {
    88           int iPivot = permuteBack[i];
    89           double sign;
    90           if (pivotRegion[i] > 0.0)
    91                sign = 1.0;
    92           else
    93                sign = -1.0;
    94           int other;
    95           if (numberInColumn[i] > 0) {
    96                int iRow = indexRow[startColumn[i]];
    97                other = permuteBack[iRow];
    98                //assert (parent_[other]!=-1);
    99           } else {
    100                other = numberRows_;
    101           }
    102           sign_[iPivot] = sign;
    103           int iParent = other;
    104           parent_[iPivot] = other;
    105           if (descendant_[iParent] >= 0) {
    106                // we have a sibling
    107                int iRight = descendant_[iParent];
    108                rightSibling_[iPivot] = iRight;
    109                leftSibling_[iRight] = iPivot;
    110           } else {
    111                rightSibling_[iPivot] = -1;
    112           }
    113           descendant_[iParent] = iPivot;
    114           leftSibling_[iPivot] = -1;
    115      }
    116      // do depth
    117      int nStack = 1;
    118      stack_[0] = descendant_[numberRows_];
    119      depth_[numberRows_] = -1; // root
    120      while (nStack) {
    121           // take off
    122           int iNext = stack_[--nStack];
    123           if (iNext >= 0) {
    124                depth_[iNext] = nStack;
    125                int iRight = rightSibling_[iNext];
    126                stack_[nStack++] = iRight;
    127                if (descendant_[iNext] >= 0)
    128                     stack_[nStack++] = descendant_[iNext];
    129           }
    130      }
    131      model_ = model;
    132      check();
     52  numberRows_ = numberRows;
     53  numberColumns_ = numberRows;
     54  parent_ = new int[numberRows_ + 1];
     55  descendant_ = new int[numberRows_ + 1];
     56  pivot_ = new int[numberRows_ + 1];
     57  rightSibling_ = new int[numberRows_ + 1];
     58  leftSibling_ = new int[numberRows_ + 1];
     59  sign_ = new double[numberRows_ + 1];
     60  stack_ = new int[numberRows_ + 1];
     61  stack2_ = new int[numberRows_ + 1];
     62  depth_ = new int[numberRows_ + 1];
     63  mark_ = new char[numberRows_ + 1];
     64  permute_ = new int[numberRows_ + 1];
     65  permuteBack_ = new int[numberRows_ + 1];
     66  int i;
     67  for (i = 0; i < numberRows_ + 1; i++) {
     68    parent_[i] = -1;
     69    descendant_[i] = -1;
     70    pivot_[i] = -1;
     71    rightSibling_[i] = -1;
     72    leftSibling_[i] = -1;
     73    sign_[i] = -1.0;
     74    stack_[i] = -1;
     75    permute_[i] = i;
     76    permuteBack_[i] = i;
     77    stack2_[i] = -1;
     78    depth_[i] = -1;
     79    mark_[i] = 0;
     80  }
     81  mark_[numberRows_] = 1;
     82  // pivotColumnBack gives order of pivoting into basis
     83  // so pivotColumnback[0] is first slack in basis and
     84  // it pivots on row permuteBack[0]
     85  // a known root is given by permuteBack[numberRows_-1]
     86  for (i = 0; i < numberRows_; i++) {
     87    int iPivot = permuteBack[i];
     88    double sign;
     89    if (pivotRegion[i] > 0.0)
     90      sign = 1.0;
     91    else
     92      sign = -1.0;
     93    int other;
     94    if (numberInColumn[i] > 0) {
     95      int iRow = indexRow[startColumn[i]];
     96      other = permuteBack[iRow];
     97      //assert (parent_[other]!=-1);
     98    } else {
     99      other = numberRows_;
     100    }
     101    sign_[iPivot] = sign;
     102    int iParent = other;
     103    parent_[iPivot] = other;
     104    if (descendant_[iParent] >= 0) {
     105      // we have a sibling
     106      int iRight = descendant_[iParent];
     107      rightSibling_[iPivot] = iRight;
     108      leftSibling_[iRight] = iPivot;
     109    } else {
     110      rightSibling_[iPivot] = -1;
     111    }
     112    descendant_[iParent] = iPivot;
     113    leftSibling_[iPivot] = -1;
     114  }
     115  // do depth
     116  int nStack = 1;
     117  stack_[0] = descendant_[numberRows_];
     118  depth_[numberRows_] = -1; // root
     119  while (nStack) {
     120    // take off
     121    int iNext = stack_[--nStack];
     122    if (iNext >= 0) {
     123      depth_[iNext] = nStack;
     124      int iRight = rightSibling_[iNext];
     125      stack_[nStack++] = iRight;
     126      if (descendant_[iNext] >= 0)
     127        stack_[nStack++] = descendant_[iNext];
     128    }
     129  }
     130  model_ = model;
     131  check();
    133132}
    134133
     
    136135// Copy constructor
    137136//-------------------------------------------------------------------
    138 ClpNetworkBasis::ClpNetworkBasis (const ClpNetworkBasis & rhs)
     137ClpNetworkBasis::ClpNetworkBasis(const ClpNetworkBasis &rhs)
    139138{
    140139#ifndef COIN_FAST_CODE
    141      slackValue_ = rhs.slackValue_;
     140  slackValue_ = rhs.slackValue_;
    142141#endif
    143      numberRows_ = rhs.numberRows_;
    144      numberColumns_ = rhs.numberColumns_;
    145      if (rhs.parent_) {
    146           parent_ = new int [numberRows_+1];
    147           CoinMemcpyN(rhs.parent_, (numberRows_ + 1), parent_);
    148      } else {
    149           parent_ = NULL;
    150      }
    151      if (rhs.descendant_) {
    152           descendant_ = new int [numberRows_+1];
    153           CoinMemcpyN(rhs.descendant_, (numberRows_ + 1), descendant_);
    154      } else {
    155           descendant_ = NULL;
    156      }
    157      if (rhs.pivot_) {
    158           pivot_ = new int [numberRows_+1];
    159           CoinMemcpyN(rhs.pivot_, (numberRows_ + 1), pivot_);
    160      } else {
    161           pivot_ = NULL;
    162      }
    163      if (rhs.rightSibling_) {
    164           rightSibling_ = new int [numberRows_+1];
    165           CoinMemcpyN(rhs.rightSibling_, (numberRows_ + 1), rightSibling_);
    166      } else {
    167           rightSibling_ = NULL;
    168      }
    169      if (rhs.leftSibling_) {
    170           leftSibling_ = new int [numberRows_+1];
    171           CoinMemcpyN(rhs.leftSibling_, (numberRows_ + 1), leftSibling_);
    172      } else {
    173           leftSibling_ = NULL;
    174      }
    175      if (rhs.sign_) {
    176           sign_ = new double [numberRows_+1];
    177           CoinMemcpyN(rhs.sign_, (numberRows_ + 1), sign_);
    178      } else {
    179           sign_ = NULL;
    180      }
    181      if (rhs.stack_) {
    182           stack_ = new int [numberRows_+1];
    183           CoinMemcpyN(rhs.stack_, (numberRows_ + 1), stack_);
    184      } else {
    185           stack_ = NULL;
    186      }
    187      if (rhs.permute_) {
    188           permute_ = new int [numberRows_+1];
    189           CoinMemcpyN(rhs.permute_, (numberRows_ + 1), permute_);
    190      } else {
    191           permute_ = NULL;
    192      }
    193      if (rhs.permuteBack_) {
    194           permuteBack_ = new int [numberRows_+1];
    195           CoinMemcpyN(rhs.permuteBack_, (numberRows_ + 1), permuteBack_);
    196      } else {
    197           permuteBack_ = NULL;
    198      }
    199      if (rhs.stack2_) {
    200           stack2_ = new int [numberRows_+1];
    201           CoinMemcpyN(rhs.stack2_, (numberRows_ + 1), stack2_);
    202      } else {
    203           stack2_ = NULL;
    204      }
    205      if (rhs.depth_) {
    206           depth_ = new int [numberRows_+1];
    207           CoinMemcpyN(rhs.depth_, (numberRows_ + 1), depth_);
    208      } else {
    209           depth_ = NULL;
    210      }
    211      if (rhs.mark_) {
    212           mark_ = new char [numberRows_+1];
    213           CoinMemcpyN(rhs.mark_, (numberRows_ + 1), mark_);
    214      } else {
    215           mark_ = NULL;
    216      }
    217      model_ = rhs.model_;
     142  numberRows_ = rhs.numberRows_;
     143  numberColumns_ = rhs.numberColumns_;
     144  if (rhs.parent_) {
     145    parent_ = new int[numberRows_ + 1];
     146    CoinMemcpyN(rhs.parent_, (numberRows_ + 1), parent_);
     147  } else {
     148    parent_ = NULL;
     149  }
     150  if (rhs.descendant_) {
     151    descendant_ = new int[numberRows_ + 1];
     152    CoinMemcpyN(rhs.descendant_, (numberRows_ + 1), descendant_);
     153  } else {
     154    descendant_ = NULL;
     155  }
     156  if (rhs.pivot_) {
     157    pivot_ = new int[numberRows_ + 1];
     158    CoinMemcpyN(rhs.pivot_, (numberRows_ + 1), pivot_);
     159  } else {
     160    pivot_ = NULL;
     161  }
     162  if (rhs.rightSibling_) {
     163    rightSibling_ = new int[numberRows_ + 1];
     164    CoinMemcpyN(rhs.rightSibling_, (numberRows_ + 1), rightSibling_);
     165  } else {
     166    rightSibling_ = NULL;
     167  }
     168  if (rhs.leftSibling_) {
     169    leftSibling_ = new int[numberRows_ + 1];
     170    CoinMemcpyN(rhs.leftSibling_, (numberRows_ + 1), leftSibling_);
     171  } else {
     172    leftSibling_ = NULL;
     173  }
     174  if (rhs.sign_) {
     175    sign_ = new double[numberRows_ + 1];
     176    CoinMemcpyN(rhs.sign_, (numberRows_ + 1), sign_);
     177  } else {
     178    sign_ = NULL;
     179  }
     180  if (rhs.stack_) {
     181    stack_ = new int[numberRows_ + 1];
     182    CoinMemcpyN(rhs.stack_, (numberRows_ + 1), stack_);
     183  } else {
     184    stack_ = NULL;
     185  }
     186  if (rhs.permute_) {
     187    permute_ = new int[numberRows_ + 1];
     188    CoinMemcpyN(rhs.permute_, (numberRows_ + 1), permute_);
     189  } else {
     190    permute_ = NULL;
     191  }
     192  if (rhs.permuteBack_) {
     193    permuteBack_ = new int[numberRows_ + 1];
     194    CoinMemcpyN(rhs.permuteBack_, (numberRows_ + 1), permuteBack_);
     195  } else {
     196    permuteBack_ = NULL;
     197  }
     198  if (rhs.stack2_) {
     199    stack2_ = new int[numberRows_ + 1];
     200    CoinMemcpyN(rhs.stack2_, (numberRows_ + 1), stack2_);
     201  } else {
     202    stack2_ = NULL;
     203  }
     204  if (rhs.depth_) {
     205    depth_ = new int[numberRows_ + 1];
     206    CoinMemcpyN(rhs.depth_, (numberRows_ + 1), depth_);
     207  } else {
     208    depth_ = NULL;
     209  }
     210  if (rhs.mark_) {
     211    mark_ = new char[numberRows_ + 1];
     212    CoinMemcpyN(rhs.mark_, (numberRows_ + 1), mark_);
     213  } else {
     214    mark_ = NULL;
     215  }
     216  model_ = rhs.model_;
    218217}
    219218
     
    221220// Destructor
    222221//-------------------------------------------------------------------
    223 ClpNetworkBasis::~ClpNetworkBasis ()
     222ClpNetworkBasis::~ClpNetworkBasis()
    224223{
    225      delete [] parent_;
    226      delete [] descendant_;
    227      delete [] pivot_;
    228      delete [] rightSibling_;
    229      delete [] leftSibling_;
    230      delete [] sign_;
    231      delete [] stack_;
    232      delete [] permute_;
    233      delete [] permuteBack_;
    234      delete [] stack2_;
    235      delete [] depth_;
    236      delete [] mark_;
     224  delete[] parent_;
     225  delete[] descendant_;
     226  delete[] pivot_;
     227  delete[] rightSibling_;
     228  delete[] leftSibling_;
     229  delete[] sign_;
     230  delete[] stack_;
     231  delete[] permute_;
     232  delete[] permuteBack_;
     233  delete[] stack2_;
     234  delete[] depth_;
     235  delete[] mark_;
    237236}
    238237
     
    241240//-------------------------------------------------------------------
    242241ClpNetworkBasis &
    243 ClpNetworkBasis::operator=(const ClpNetworkBasis& rhs)
     242ClpNetworkBasis::operator=(const ClpNetworkBasis &rhs)
    244243{
    245      if (this != &rhs) {
    246           delete [] parent_;
    247           delete [] descendant_;
    248           delete [] pivot_;
    249           delete [] rightSibling_;
    250           delete [] leftSibling_;
    251           delete [] sign_;
    252           delete [] stack_;
    253           delete [] permute_;
    254           delete [] permuteBack_;
    255           delete [] stack2_;
    256           delete [] depth_;
    257           delete [] mark_;
     244  if (this != &rhs) {
     245    delete[] parent_;
     246    delete[] descendant_;
     247    delete[] pivot_;
     248    delete[] rightSibling_;
     249    delete[] leftSibling_;
     250    delete[] sign_;
     251    delete[] stack_;
     252    delete[] permute_;
     253    delete[] permuteBack_;
     254    delete[] stack2_;
     255    delete[] depth_;
     256    delete[] mark_;
    258257#ifndef COIN_FAST_CODE
    259           slackValue_ = rhs.slackValue_;
     258    slackValue_ = rhs.slackValue_;
    260259#endif
    261           numberRows_ = rhs.numberRows_;
    262           numberColumns_ = rhs.numberColumns_;
    263           if (rhs.parent_) {
    264                parent_ = new int [numberRows_+1];
    265                CoinMemcpyN(rhs.parent_, (numberRows_ + 1), parent_);
    266           } else {
    267                parent_ = NULL;
    268           }
    269           if (rhs.descendant_) {
    270                descendant_ = new int [numberRows_+1];
    271                CoinMemcpyN(rhs.descendant_, (numberRows_ + 1), descendant_);
    272           } else {
    273                descendant_ = NULL;
    274           }
    275           if (rhs.pivot_) {
    276                pivot_ = new int [numberRows_+1];
    277                CoinMemcpyN(rhs.pivot_, (numberRows_ + 1), pivot_);
    278           } else {
    279                pivot_ = NULL;
    280           }
    281           if (rhs.rightSibling_) {
    282                rightSibling_ = new int [numberRows_+1];
    283                CoinMemcpyN(rhs.rightSibling_, (numberRows_ + 1), rightSibling_);
    284           } else {
    285                rightSibling_ = NULL;
    286           }
    287           if (rhs.leftSibling_) {
    288                leftSibling_ = new int [numberRows_+1];
    289                CoinMemcpyN(rhs.leftSibling_, (numberRows_ + 1), leftSibling_);
    290           } else {
    291                leftSibling_ = NULL;
    292           }
    293           if (rhs.sign_) {
    294                sign_ = new double [numberRows_+1];
    295                CoinMemcpyN(rhs.sign_, (numberRows_ + 1), sign_);
    296           } else {
    297                sign_ = NULL;
    298           }
    299           if (rhs.stack_) {
    300                stack_ = new int [numberRows_+1];
    301                CoinMemcpyN(rhs.stack_, (numberRows_ + 1), stack_);
    302           } else {
    303                stack_ = NULL;
    304           }
    305           if (rhs.permute_) {
    306                permute_ = new int [numberRows_+1];
    307                CoinMemcpyN(rhs.permute_, (numberRows_ + 1), permute_);
    308           } else {
    309                permute_ = NULL;
    310           }
    311           if (rhs.permuteBack_) {
    312                permuteBack_ = new int [numberRows_+1];
    313                CoinMemcpyN(rhs.permuteBack_, (numberRows_ + 1), permuteBack_);
    314           } else {
    315                permuteBack_ = NULL;
    316           }
    317           if (rhs.stack2_) {
    318                stack2_ = new int [numberRows_+1];
    319                CoinMemcpyN(rhs.stack2_, (numberRows_ + 1), stack2_);
    320           } else {
    321                stack2_ = NULL;
    322           }
    323           if (rhs.depth_) {
    324                depth_ = new int [numberRows_+1];
    325                CoinMemcpyN(rhs.depth_, (numberRows_ + 1), depth_);
    326           } else {
    327                depth_ = NULL;
    328           }
    329           if (rhs.mark_) {
    330                mark_ = new char [numberRows_+1];
    331                CoinMemcpyN(rhs.mark_, (numberRows_ + 1), mark_);
    332           } else {
    333                mark_ = NULL;
    334           }
    335      }
    336      return *this;
     260    numberRows_ = rhs.numberRows_;
     261    numberColumns_ = rhs.numberColumns_;
     262    if (rhs.parent_) {
     263      parent_ = new int[numberRows_ + 1];
     264      CoinMemcpyN(rhs.parent_, (numberRows_ + 1), parent_);
     265    } else {
     266      parent_ = NULL;
     267    }
     268    if (rhs.descendant_) {
     269      descendant_ = new int[numberRows_ + 1];
     270      CoinMemcpyN(rhs.descendant_, (numberRows_ + 1), descendant_);
     271    } else {
     272      descendant_ = NULL;
     273    }
     274    if (rhs.pivot_) {
     275      pivot_ = new int[numberRows_ + 1];
     276      CoinMemcpyN(rhs.pivot_, (numberRows_ + 1), pivot_);
     277    } else {
     278      pivot_ = NULL;
     279    }
     280    if (rhs.rightSibling_) {
     281      rightSibling_ = new int[numberRows_ + 1];
     282      CoinMemcpyN(rhs.rightSibling_, (numberRows_ + 1), rightSibling_);
     283    } else {
     284      rightSibling_ = NULL;
     285    }
     286    if (rhs.leftSibling_) {
     287      leftSibling_ = new int[numberRows_ + 1];
     288      CoinMemcpyN(rhs.leftSibling_, (numberRows_ + 1), leftSibling_);
     289    } else {
     290      leftSibling_ = NULL;
     291    }
     292    if (rhs.sign_) {
     293      sign_ = new double[numberRows_ + 1];
     294      CoinMemcpyN(rhs.sign_, (numberRows_ + 1), sign_);
     295    } else {
     296      sign_ = NULL;
     297    }
     298    if (rhs.stack_) {
     299      stack_ = new int[numberRows_ + 1];
     300      CoinMemcpyN(rhs.stack_, (numberRows_ + 1), stack_);
     301    } else {
     302      stack_ = NULL;
     303    }
     304    if (rhs.permute_) {
     305      permute_ = new int[numberRows_ + 1];
     306      CoinMemcpyN(rhs.permute_, (numberRows_ + 1), permute_);
     307    } else {
     308      permute_ = NULL;
     309    }
     310    if (rhs.permuteBack_) {
     311      permuteBack_ = new int[numberRows_ + 1];
     312      CoinMemcpyN(rhs.permuteBack_, (numberRows_ + 1), permuteBack_);
     313    } else {
     314      permuteBack_ = NULL;
     315    }
     316    if (rhs.stack2_) {
     317      stack2_ = new int[numberRows_ + 1];
     318      CoinMemcpyN(rhs.stack2_, (numberRows_ + 1), stack2_);
     319    } else {
     320      stack2_ = NULL;
     321    }
     322    if (rhs.depth_) {
     323      depth_ = new int[numberRows_ + 1];
     324      CoinMemcpyN(rhs.depth_, (numberRows_ + 1), depth_);
     325    } else {
     326      depth_ = NULL;
     327    }
     328    if (rhs.mark_) {
     329      mark_ = new char[numberRows_ + 1];
     330      CoinMemcpyN(rhs.mark_, (numberRows_ + 1), mark_);
     331    } else {
     332      mark_ = NULL;
     333    }
     334  }
     335  return *this;
    337336}
    338337// checks looks okay
    339338void ClpNetworkBasis::check()
    340339{
    341      // check depth
    342      int nStack = 1;
    343      stack_[0] = descendant_[numberRows_];
    344      depth_[numberRows_] = -1; // root
    345      while (nStack) {
    346           // take off
    347           int iNext = stack_[--nStack];
    348           if (iNext >= 0) {
    349                //assert (depth_[iNext]==nStack);
    350                depth_[iNext] = nStack;
    351                int iRight = rightSibling_[iNext];
    352                stack_[nStack++] = iRight;
    353                if (descendant_[iNext] >= 0)
    354                     stack_[nStack++] = descendant_[iNext];
    355           }
    356      }
     340  // check depth
     341  int nStack = 1;
     342  stack_[0] = descendant_[numberRows_];
     343  depth_[numberRows_] = -1; // root
     344  while (nStack) {
     345    // take off
     346    int iNext = stack_[--nStack];
     347    if (iNext >= 0) {
     348      //assert (depth_[iNext]==nStack);
     349      depth_[iNext] = nStack;
     350      int iRight = rightSibling_[iNext];
     351      stack_[nStack++] = iRight;
     352      if (descendant_[iNext] >= 0)
     353        stack_[nStack++] = descendant_[iNext];
     354    }
     355  }
    357356}
    358357// prints
    359358void ClpNetworkBasis::print()
    360359{
    361      int i;
    362      printf("       parent descendant     left    right   sign    depth\n");
    363      for (i = 0; i < numberRows_ + 1; i++)
    364           printf("%4d  %7d   %8d  %7d  %7d  %5g  %7d\n",
    365                  i, parent_[i], descendant_[i], leftSibling_[i], rightSibling_[i],
    366                  sign_[i], depth_[i]);
     360  int i;
     361  printf("       parent descendant     left    right   sign    depth\n");
     362  for (i = 0; i < numberRows_ + 1; i++)
     363    printf("%4d  %7d   %8d  %7d  %7d  %5g  %7d\n",
     364      i, parent_[i], descendant_[i], leftSibling_[i], rightSibling_[i],
     365      sign_[i], depth_[i]);
    367366}
    368367/* Replaces one Column to basis,
    369368   returns 0=OK
    370369*/
    371 int
    372 ClpNetworkBasis::replaceColumn ( CoinIndexedVector * regionSparse,
    373                                  int pivotRow)
     370int ClpNetworkBasis::replaceColumn(CoinIndexedVector *regionSparse,
     371  int pivotRow)
    374372{
    375      // When things have settled down then redo this to make more elegant
    376      // I am sure lots of loops can be combined
    377      // regionSparse is empty
    378      assert (!regionSparse->getNumElements());
    379      model_->unpack(regionSparse, model_->sequenceIn());
    380      // arc given by pivotRow is leaving basis
    381      //int kParent = parent_[pivotRow];
    382      // arc coming in has these two nodes
    383      int * indices = regionSparse->getIndices();
    384      int iRow0 = indices[0];
    385      int iRow1;
    386      if (regionSparse->getNumElements() == 2)
    387           iRow1 = indices[1];
    388      else
    389           iRow1 = numberRows_;
    390      double sign = -regionSparse->denseVector()[iRow0];
    391      regionSparse->clear();
    392      // and outgoing
    393      model_->unpack(regionSparse, model_->pivotVariable()[pivotRow]);
    394      int jRow0 = indices[0];
    395      int jRow1;
    396      if (regionSparse->getNumElements() == 2)
    397           jRow1 = indices[1];
    398      else
    399           jRow1 = numberRows_;
    400      regionSparse->clear();
    401      // get correct pivotRow
    402      //#define FULL_DEBUG
     373  // When things have settled down then redo this to make more elegant
     374  // I am sure lots of loops can be combined
     375  // regionSparse is empty
     376  assert(!regionSparse->getNumElements());
     377  model_->unpack(regionSparse, model_->sequenceIn());
     378  // arc given by pivotRow is leaving basis
     379  //int kParent = parent_[pivotRow];
     380  // arc coming in has these two nodes
     381  int *indices = regionSparse->getIndices();
     382  int iRow0 = indices[0];
     383  int iRow1;
     384  if (regionSparse->getNumElements() == 2)
     385    iRow1 = indices[1];
     386  else
     387    iRow1 = numberRows_;
     388  double sign = -regionSparse->denseVector()[iRow0];
     389  regionSparse->clear();
     390  // and outgoing
     391  model_->unpack(regionSparse, model_->pivotVariable()[pivotRow]);
     392  int jRow0 = indices[0];
     393  int jRow1;
     394  if (regionSparse->getNumElements() == 2)
     395    jRow1 = indices[1];
     396  else
     397    jRow1 = numberRows_;
     398  regionSparse->clear();
     399  // get correct pivotRow
     400  //#define FULL_DEBUG
    403401#ifdef FULL_DEBUG
    404      printf ("irow %d %d, jrow %d %d\n",
    405              iRow0, iRow1, jRow0, jRow1);
     402  printf("irow %d %d, jrow %d %d\n",
     403    iRow0, iRow1, jRow0, jRow1);
    406404#endif
    407      if (parent_[jRow0] == jRow1) {
    408           int newPivot = jRow0;
    409           if (newPivot != pivotRow) {
     405  if (parent_[jRow0] == jRow1) {
     406    int newPivot = jRow0;
     407    if (newPivot != pivotRow) {
    410408#ifdef FULL_DEBUG
    411                printf("pivot row of %d permuted to %d\n", pivotRow, newPivot);
     409      printf("pivot row of %d permuted to %d\n", pivotRow, newPivot);
    412410#endif
    413                pivotRow = newPivot;
    414           }
    415      } else {
    416           //assert (parent_[jRow1]==jRow0);
    417           int newPivot = jRow1;
    418           if (newPivot != pivotRow) {
     411      pivotRow = newPivot;
     412    }
     413  } else {
     414    //assert (parent_[jRow1]==jRow0);
     415    int newPivot = jRow1;
     416    if (newPivot != pivotRow) {
    419417#ifdef FULL_DEBUG
    420                printf("pivot row of %d permuted to %d\n", pivotRow, newPivot);
     418      printf("pivot row of %d permuted to %d\n", pivotRow, newPivot);
    421419#endif
    422                pivotRow = newPivot;
    423           }
    424      }
    425      bool extraPrint = (model_->numberIterations() > -3) &&
    426                        (model_->logLevel() > 10);
    427      if (extraPrint)
    428           print();
     420      pivotRow = newPivot;
     421    }
     422  }
     423  bool extraPrint = (model_->numberIterations() > -3) && (model_->logLevel() > 10);
     424  if (extraPrint)
     425    print();
    429426#ifdef FULL_DEBUG
    430      printf("In %d (region= %g, stored %g) %d (%g) pivoting on %d (%g)\n",
    431             iRow1, sign, sign_[iRow1], iRow0, sign_[iRow0] , pivotRow, sign_[pivotRow]);
     427  printf("In %d (region= %g, stored %g) %d (%g) pivoting on %d (%g)\n",
     428    iRow1, sign, sign_[iRow1], iRow0, sign_[iRow0], pivotRow, sign_[pivotRow]);
    432429#endif
    433      // see which path outgoing pivot is on
    434      int kRow = -1;
    435      int jRow = iRow1;
    436      while (jRow != numberRows_) {
    437           if (jRow == pivotRow) {
    438                kRow = iRow1;
    439                break;
    440           } else {
    441                jRow = parent_[jRow];
    442           }
    443      }
    444      if (kRow < 0) {
    445           jRow = iRow0;
    446           while (jRow != numberRows_) {
    447                if (jRow == pivotRow) {
    448                     kRow = iRow0;
    449                     break;
    450                } else {
    451                     jRow = parent_[jRow];
    452                }
    453           }
    454      }
    455      //assert (kRow>=0);
    456      if (iRow0 == kRow) {
    457           iRow0 = iRow1;
    458           iRow1 = kRow;
    459           sign = -sign;
    460      }
    461      // pivot row is on path from iRow1 back to root
    462      // get stack of nodes to change
    463      // Also get precursors for cleaning order
    464      int nStack = 1;
    465      stack_[0] = iRow0;
    466      while (kRow != pivotRow) {
    467           stack_[nStack++] = kRow;
    468           if (sign * sign_[kRow] < 0.0) {
    469                sign_[kRow] = -sign_[kRow];
    470           } else {
    471                sign = -sign;
    472           }
    473           kRow = parent_[kRow];
    474           //sign *= sign_[kRow];
    475      }
    476      stack_[nStack++] = pivotRow;
    477      if (sign * sign_[pivotRow] < 0.0) {
    478           sign_[pivotRow] = -sign_[pivotRow];
    479      } else {
    480           sign = -sign;
    481      }
    482      int iParent = parent_[pivotRow];
    483      while (nStack > 1) {
    484           int iLeft;
    485           int iRight;
    486           kRow = stack_[--nStack];
    487           int newParent = stack_[nStack-1];
     430  // see which path outgoing pivot is on
     431  int kRow = -1;
     432  int jRow = iRow1;
     433  while (jRow != numberRows_) {
     434    if (jRow == pivotRow) {
     435      kRow = iRow1;
     436      break;
     437    } else {
     438      jRow = parent_[jRow];
     439    }
     440  }
     441  if (kRow < 0) {
     442    jRow = iRow0;
     443    while (jRow != numberRows_) {
     444      if (jRow == pivotRow) {
     445        kRow = iRow0;
     446        break;
     447      } else {
     448        jRow = parent_[jRow];
     449      }
     450    }
     451  }
     452  //assert (kRow>=0);
     453  if (iRow0 == kRow) {
     454    iRow0 = iRow1;
     455    iRow1 = kRow;
     456    sign = -sign;
     457  }
     458  // pivot row is on path from iRow1 back to root
     459  // get stack of nodes to change
     460  // Also get precursors for cleaning order
     461  int nStack = 1;
     462  stack_[0] = iRow0;
     463  while (kRow != pivotRow) {
     464    stack_[nStack++] = kRow;
     465    if (sign * sign_[kRow] < 0.0) {
     466      sign_[kRow] = -sign_[kRow];
     467    } else {
     468      sign = -sign;
     469    }
     470    kRow = parent_[kRow];
     471    //sign *= sign_[kRow];
     472  }
     473  stack_[nStack++] = pivotRow;
     474  if (sign * sign_[pivotRow] < 0.0) {
     475    sign_[pivotRow] = -sign_[pivotRow];
     476  } else {
     477    sign = -sign;
     478  }
     479  int iParent = parent_[pivotRow];
     480  while (nStack > 1) {
     481    int iLeft;
     482    int iRight;
     483    kRow = stack_[--nStack];
     484    int newParent = stack_[nStack - 1];
    488485#ifdef FULL_DEBUG
    489           printf("row %d, old parent %d, new parent %d, pivotrow %d\n", kRow,
    490                  iParent, newParent, pivotRow);
     486    printf("row %d, old parent %d, new parent %d, pivotrow %d\n", kRow,
     487      iParent, newParent, pivotRow);
    491488#endif
    492           int i1 = permuteBack_[pivotRow];
    493           int i2 = permuteBack_[kRow];
    494           permuteBack_[pivotRow] = i2;
    495           permuteBack_[kRow] = i1;
    496           // do Btran permutation
    497           permute_[i1] = kRow;
    498           permute_[i2] = pivotRow;
    499           pivotRow = kRow;
    500           // Take out of old parent
    501           iLeft = leftSibling_[kRow];
    502           iRight = rightSibling_[kRow];
    503           if (iLeft < 0) {
    504                // take out of tree
    505                if (iRight >= 0) {
    506                     leftSibling_[iRight] = iLeft;
    507                     descendant_[iParent] = iRight;
    508                } else {
     489    int i1 = permuteBack_[pivotRow];
     490    int i2 = permuteBack_[kRow];
     491    permuteBack_[pivotRow] = i2;
     492    permuteBack_[kRow] = i1;
     493    // do Btran permutation
     494    permute_[i1] = kRow;
     495    permute_[i2] = pivotRow;
     496    pivotRow = kRow;
     497    // Take out of old parent
     498    iLeft = leftSibling_[kRow];
     499    iRight = rightSibling_[kRow];
     500    if (iLeft < 0) {
     501      // take out of tree
     502      if (iRight >= 0) {
     503        leftSibling_[iRight] = iLeft;
     504        descendant_[iParent] = iRight;
     505      } else {
    509506#ifdef FULL_DEBUG
    510                     printf("Saying %d (old parent of %d) has no descendants\n",
    511                            iParent, kRow);
     507        printf("Saying %d (old parent of %d) has no descendants\n",
     508          iParent, kRow);
    512509#endif
    513                     descendant_[iParent] = -1;
    514                }
    515           } else {
    516                // take out of tree
    517                rightSibling_[iLeft] = iRight;
    518                if (iRight >= 0)
    519                     leftSibling_[iRight] = iLeft;
    520           }
    521           leftSibling_[kRow] = -1;
    522           rightSibling_[kRow] = -1;
     510        descendant_[iParent] = -1;
     511      }
     512    } else {
     513      // take out of tree
     514      rightSibling_[iLeft] = iRight;
     515      if (iRight >= 0)
     516        leftSibling_[iRight] = iLeft;
     517    }
     518    leftSibling_[kRow] = -1;
     519    rightSibling_[kRow] = -1;
    523520
    524           // now insert new one
    525           // make this descendant of that
    526           if (descendant_[newParent] >= 0) {
    527                // we will have a sibling
    528                int jRight = descendant_[newParent];
    529                rightSibling_[kRow] = jRight;
    530                leftSibling_[jRight] = kRow;
    531           } else {
    532                rightSibling_[kRow] = -1;
    533           }
    534           descendant_[newParent] = kRow;
    535           leftSibling_[kRow] = -1;
    536           parent_[kRow] = newParent;
     521    // now insert new one
     522    // make this descendant of that
     523    if (descendant_[newParent] >= 0) {
     524      // we will have a sibling
     525      int jRight = descendant_[newParent];
     526      rightSibling_[kRow] = jRight;
     527      leftSibling_[jRight] = kRow;
     528    } else {
     529      rightSibling_[kRow] = -1;
     530    }
     531    descendant_[newParent] = kRow;
     532    leftSibling_[kRow] = -1;
     533    parent_[kRow] = newParent;
    537534
    538           iParent = kRow;
    539      }
    540      // now redo all depths from stack_[1]
    541      // This must be possible to combine - but later
    542      {
    543           int iPivot = stack_[1];
    544           int iDepth = depth_[parent_[iPivot]]; //depth of parent
    545           iDepth ++; //correct for this one
    546           int nStack = 1;
    547           stack_[0] = iPivot;
    548           while (nStack) {
    549                // take off
    550                int iNext = stack_[--nStack];
    551                if (iNext >= 0) {
    552                     // add stack level
    553                     depth_[iNext] = nStack + iDepth;
    554                     stack_[nStack++] = rightSibling_[iNext];
    555                     if (descendant_[iNext] >= 0)
    556                          stack_[nStack++] = descendant_[iNext];
    557                }
    558           }
    559      }
    560      if (extraPrint)
    561           print();
    562      //check();
    563      return 0;
     535    iParent = kRow;
     536  }
     537  // now redo all depths from stack_[1]
     538  // This must be possible to combine - but later
     539  {
     540    int iPivot = stack_[1];
     541    int iDepth = depth_[parent_[iPivot]]; //depth of parent
     542    iDepth++; //correct for this one
     543    int nStack = 1;
     544    stack_[0] = iPivot;
     545    while (nStack) {
     546      // take off
     547      int iNext = stack_[--nStack];
     548      if (iNext >= 0) {
     549        // add stack level
     550        depth_[iNext] = nStack + iDepth;
     551        stack_[nStack++] = rightSibling_[iNext];
     552        if (descendant_[iNext] >= 0)
     553          stack_[nStack++] = descendant_[iNext];
     554      }
     555    }
     556  }
     557  if (extraPrint)
     558    print();
     559  //check();
     560  return 0;
    564561}
    565562
    566563/* Updates one column (FTRAN) from region2 */
    567564double
    568 ClpNetworkBasis::updateColumn (  CoinIndexedVector * regionSparse,
    569                                  CoinIndexedVector * regionSparse2,
    570                                  int pivotRow)
     565ClpNetworkBasis::updateColumn(CoinIndexedVector *regionSparse,
     566  CoinIndexedVector *regionSparse2,
     567  int pivotRow)
    571568{
    572      regionSparse->clear (  );
    573      double *region = regionSparse->denseVector (  );
    574      double *region2 = regionSparse2->denseVector (  );
    575      int *regionIndex2 = regionSparse2->getIndices (  );
    576      int numberNonZero = regionSparse2->getNumElements (  );
    577      int *regionIndex = regionSparse->getIndices (  );
    578      int i;
    579      bool doTwo = (numberNonZero == 2);
    580      int i0 = -1;
    581      int i1 = -1;
    582      if (doTwo) {
    583           i0 = regionIndex2[0];
    584           i1 = regionIndex2[1];
    585      }
    586      double returnValue = 0.0;
    587      bool packed = regionSparse2->packedMode();
    588      if (packed) {
    589           if (doTwo && region2[0]*region2[1] < 0.0) {
    590                region[i0] = region2[0];
    591                region2[0] = 0.0;
    592                region[i1] = region2[1];
    593                region2[1] = 0.0;
    594                int iDepth0 = depth_[i0];
    595                int iDepth1 = depth_[i1];
    596                if (iDepth1 > iDepth0) {
    597                     int temp = i0;
    598                     i0 = i1;
    599                     i1 = temp;
    600                     temp = iDepth0;
    601                     iDepth0 = iDepth1;
    602                     iDepth1 = temp;
    603                }
    604                numberNonZero = 0;
    605                if (pivotRow < 0) {
    606                     while (iDepth0 > iDepth1) {
    607                          double pivotValue = region[i0];
    608                          // put back now ?
    609                          int iBack = permuteBack_[i0];
    610                          region2[numberNonZero] = pivotValue * sign_[i0];
    611                          regionIndex2[numberNonZero++] = iBack;
    612                          int otherRow = parent_[i0];
    613                          region[i0] = 0.0;
    614                          region[otherRow] += pivotValue;
    615                          iDepth0--;
    616                          i0 = otherRow;
    617                     }
    618                     while (i0 != i1) {
    619                          double pivotValue = region[i0];
    620                          // put back now ?
    621                          int iBack = permuteBack_[i0];
    622                          region2[numberNonZero] = pivotValue * sign_[i0];
    623                          regionIndex2[numberNonZero++] = iBack;
    624                          int otherRow = parent_[i0];
    625                          region[i0] = 0.0;
    626                          region[otherRow] += pivotValue;
    627                          i0 = otherRow;
    628                          double pivotValue1 = region[i1];
    629                          // put back now ?
    630                          int iBack1 = permuteBack_[i1];
    631                          region2[numberNonZero] = pivotValue1 * sign_[i1];
    632                          regionIndex2[numberNonZero++] = iBack1;
    633                          int otherRow1 = parent_[i1];
    634                          region[i1] = 0.0;
    635                          region[otherRow1] += pivotValue1;
    636                          i1 = otherRow1;
    637                     }
    638                } else {
    639                     while (iDepth0 > iDepth1) {
    640                          double pivotValue = region[i0];
    641                          // put back now ?
    642                          int iBack = permuteBack_[i0];
    643                          double value = pivotValue * sign_[i0];
    644                          region2[numberNonZero] = value;
    645                          regionIndex2[numberNonZero++] = iBack;
    646                          if (iBack == pivotRow)
    647                               returnValue = value;
    648                          int otherRow = parent_[i0];
    649                          region[i0] = 0.0;
    650                          region[otherRow] += pivotValue;
    651                          iDepth0--;
    652                          i0 = otherRow;
    653                     }
    654                     while (i0 != i1) {
    655                          double pivotValue = region[i0];
    656                          // put back now ?
    657                          int iBack = permuteBack_[i0];
    658                          double value = pivotValue * sign_[i0];
    659                          region2[numberNonZero] = value;
    660                          regionIndex2[numberNonZero++] = iBack;
    661                          if (iBack == pivotRow)
    662                               returnValue = value;
    663                          int otherRow = parent_[i0];
    664                          region[i0] = 0.0;
    665                          region[otherRow] += pivotValue;
    666                          i0 = otherRow;
    667                          double pivotValue1 = region[i1];
    668                          // put back now ?
    669                          int iBack1 = permuteBack_[i1];
    670                          value = pivotValue1 * sign_[i1];
    671                          region2[numberNonZero] = value;
    672                          regionIndex2[numberNonZero++] = iBack1;
    673                          if (iBack1 == pivotRow)
    674                               returnValue = value;
    675                          int otherRow1 = parent_[i1];
    676                          region[i1] = 0.0;
    677                          region[otherRow1] += pivotValue1;
    678                          i1 = otherRow1;
    679                     }
    680                }
    681           } else {
    682                // set up linked lists at each depth
    683                // stack2 is start, stack is next
    684                int greatestDepth = -1;
    685                //mark_[numberRows_]=1;
    686                for (i = 0; i < numberNonZero; i++) {
    687                     int j = regionIndex2[i];
    688                     double value = region2[i];
    689                     region2[i] = 0.0;
    690                     region[j] = value;
    691                     regionIndex[i] = j;
    692                     int iDepth = depth_[j];
    693                     if (iDepth > greatestDepth)
    694                          greatestDepth = iDepth;
    695                     // and back until marked
    696                     while (!mark_[j]) {
    697                          int iNext = stack2_[iDepth];
    698                          stack2_[iDepth] = j;
    699                          stack_[j] = iNext;
    700                          mark_[j] = 1;
    701                          iDepth--;
    702                          j = parent_[j];
    703                     }
    704                }
    705                numberNonZero = 0;
    706                if (pivotRow < 0) {
    707                     for (; greatestDepth >= 0; greatestDepth--) {
    708                          int iPivot = stack2_[greatestDepth];
    709                          stack2_[greatestDepth] = -1;
    710                          while (iPivot >= 0) {
    711                               mark_[iPivot] = 0;
    712                               double pivotValue = region[iPivot];
    713                               if (pivotValue) {
    714                                    // put back now ?
    715                                    int iBack = permuteBack_[iPivot];
    716                                    region2[numberNonZero] = pivotValue * sign_[iPivot];
    717                                    regionIndex2[numberNonZero++] = iBack;
    718                                    int otherRow = parent_[iPivot];
    719                                    region[iPivot] = 0.0;
    720                                    region[otherRow] += pivotValue;
    721                               }
    722                               iPivot = stack_[iPivot];
    723                          }
    724                     }
    725                } else {
    726                     for (; greatestDepth >= 0; greatestDepth--) {
    727                          int iPivot = stack2_[greatestDepth];
    728                          stack2_[greatestDepth] = -1;
    729                          while (iPivot >= 0) {
    730                               mark_[iPivot] = 0;
    731                               double pivotValue = region[iPivot];
    732                               if (pivotValue) {
    733                                    // put back now ?
    734                                    int iBack = permuteBack_[iPivot];
    735                                    double value = pivotValue * sign_[iPivot];
    736                                    region2[numberNonZero] = value;
    737                                    regionIndex2[numberNonZero++] = iBack;
    738                                    if (iBack == pivotRow)
    739                                         returnValue = value;
    740                                    int otherRow = parent_[iPivot];
    741                                    region[iPivot] = 0.0;
    742                                    region[otherRow] += pivotValue;
    743                               }
    744                               iPivot = stack_[iPivot];
    745                          }
    746                     }
    747                }
     569  regionSparse->clear();
     570  double *region = regionSparse->denseVector();
     571  double *region2 = regionSparse2->denseVector();
     572  int *regionIndex2 = regionSparse2->getIndices();
     573  int numberNonZero = regionSparse2->getNumElements();
     574  int *regionIndex = regionSparse->getIndices();
     575  int i;
     576  bool doTwo = (numberNonZero == 2);
     577  int i0 = -1;
     578  int i1 = -1;
     579  if (doTwo) {
     580    i0 = regionIndex2[0];
     581    i1 = regionIndex2[1];
     582  }
     583  double returnValue = 0.0;
     584  bool packed = regionSparse2->packedMode();
     585  if (packed) {
     586    if (doTwo && region2[0] * region2[1] < 0.0) {
     587      region[i0] = region2[0];
     588      region2[0] = 0.0;
     589      region[i1] = region2[1];
     590      region2[1] = 0.0;
     591      int iDepth0 = depth_[i0];
     592      int iDepth1 = depth_[i1];
     593      if (iDepth1 > iDepth0) {
     594        int temp = i0;
     595        i0 = i1;
     596        i1 = temp;
     597        temp = iDepth0;
     598        iDepth0 = iDepth1;
     599        iDepth1 = temp;
     600      }
     601      numberNonZero = 0;
     602      if (pivotRow < 0) {
     603        while (iDepth0 > iDepth1) {
     604          double pivotValue = region[i0];
     605          // put back now ?
     606          int iBack = permuteBack_[i0];
     607          region2[numberNonZero] = pivotValue * sign_[i0];
     608          regionIndex2[numberNonZero++] = iBack;
     609          int otherRow = parent_[i0];
     610          region[i0] = 0.0;
     611          region[otherRow] += pivotValue;
     612          iDepth0--;
     613          i0 = otherRow;
     614        }
     615        while (i0 != i1) {
     616          double pivotValue = region[i0];
     617          // put back now ?
     618          int iBack = permuteBack_[i0];
     619          region2[numberNonZero] = pivotValue * sign_[i0];
     620          regionIndex2[numberNonZero++] = iBack;
     621          int otherRow = parent_[i0];
     622          region[i0] = 0.0;
     623          region[otherRow] += pivotValue;
     624          i0 = otherRow;
     625          double pivotValue1 = region[i1];
     626          // put back now ?
     627          int iBack1 = permuteBack_[i1];
     628          region2[numberNonZero] = pivotValue1 * sign_[i1];
     629          regionIndex2[numberNonZero++] = iBack1;
     630          int otherRow1 = parent_[i1];
     631          region[i1] = 0.0;
     632          region[otherRow1] += pivotValue1;
     633          i1 = otherRow1;
     634        }
     635      } else {
     636        while (iDepth0 > iDepth1) {
     637          double pivotValue = region[i0];
     638          // put back now ?
     639          int iBack = permuteBack_[i0];
     640          double value = pivotValue * sign_[i0];
     641          region2[numberNonZero] = value;
     642          regionIndex2[numberNonZero++] = iBack;
     643          if (iBack == pivotRow)
     644            returnValue = value;
     645          int otherRow = parent_[i0];
     646          region[i0] = 0.0;
     647          region[otherRow] += pivotValue;
     648          iDepth0--;
     649          i0 = otherRow;
     650        }
     651        while (i0 != i1) {
     652          double pivotValue = region[i0];
     653          // put back now ?
     654          int iBack = permuteBack_[i0];
     655          double value = pivotValue * sign_[i0];
     656          region2[numberNonZero] = value;
     657          regionIndex2[numberNonZero++] = iBack;
     658          if (iBack == pivotRow)
     659            returnValue = value;
     660          int otherRow = parent_[i0];
     661          region[i0] = 0.0;
     662          region[otherRow] += pivotValue;
     663          i0 = otherRow;
     664          double pivotValue1 = region[i1];
     665          // put back now ?
     666          int iBack1 = permuteBack_[i1];
     667          value = pivotValue1 * sign_[i1];
     668          region2[numberNonZero] = value;
     669          regionIndex2[numberNonZero++] = iBack1;
     670          if (iBack1 == pivotRow)
     671            returnValue = value;
     672          int otherRow1 = parent_[i1];
     673          region[i1] = 0.0;
     674          region[otherRow1] += pivotValue1;
     675          i1 = otherRow1;
     676        }
     677      }
     678    } else {
     679      // set up linked lists at each depth
     680      // stack2 is start, stack is next
     681      int greatestDepth = -1;
     682      //mark_[numberRows_]=1;
     683      for (i = 0; i < numberNonZero; i++) {
     684        int j = regionIndex2[i];
     685        double value = region2[i];
     686        region2[i] = 0.0;
     687        region[j] = value;
     688        regionIndex[i] = j;
     689        int iDepth = depth_[j];
     690        if (iDepth > greatestDepth)
     691          greatestDepth = iDepth;
     692        // and back until marked
     693        while (!mark_[j]) {
     694          int iNext = stack2_[iDepth];
     695          stack2_[iDepth] = j;
     696          stack_[j] = iNext;
     697          mark_[j] = 1;
     698          iDepth--;
     699          j = parent_[j];
     700        }
     701      }
     702      numberNonZero = 0;
     703      if (pivotRow < 0) {
     704        for (; greatestDepth >= 0; greatestDepth--) {
     705          int iPivot = stack2_[greatestDepth];
     706          stack2_[greatestDepth] = -1;
     707          while (iPivot >= 0) {
     708            mark_[iPivot] = 0;
     709            double pivotValue = region[iPivot];
     710            if (pivotValue) {
     711              // put back now ?
     712              int iBack = permuteBack_[iPivot];
     713              region2[numberNonZero] = pivotValue * sign_[iPivot];
     714              regionIndex2[numberNonZero++] = iBack;
     715              int otherRow = parent_[iPivot];
     716              region[iPivot] = 0.0;
     717              region[otherRow] += pivotValue;
     718            }
     719            iPivot = stack_[iPivot];
    748720          }
    749      } else {
    750           if (doTwo && region2[i0]*region2[i1] < 0.0) {
    751                // If just +- 1 then could go backwards on depth until join
    752                region[i0] = region2[i0];
    753                region2[i0] = 0.0;
    754                region[i1] = region2[i1];
    755                region2[i1] = 0.0;
    756                int iDepth0 = depth_[i0];
    757                int iDepth1 = depth_[i1];
    758                if (iDepth1 > iDepth0) {
    759                     int temp = i0;
    760                     i0 = i1;
    761                     i1 = temp;
    762                     temp = iDepth0;
    763                     iDepth0 = iDepth1;
    764                     iDepth1 = temp;
    765                }
    766                numberNonZero = 0;
    767                while (iDepth0 > iDepth1) {
    768                     double pivotValue = region[i0];
    769                     // put back now ?
    770                     int iBack = permuteBack_[i0];
    771                     regionIndex2[numberNonZero++] = iBack;
    772                     int otherRow = parent_[i0];
    773                     region2[iBack] = pivotValue * sign_[i0];
    774                     region[i0] = 0.0;
    775                     region[otherRow] += pivotValue;
    776                     iDepth0--;
    777                     i0 = otherRow;
    778                }
    779                while (i0 != i1) {
    780                     double pivotValue = region[i0];
    781                     // put back now ?
    782                     int iBack = permuteBack_[i0];
    783                     regionIndex2[numberNonZero++] = iBack;
    784                     int otherRow = parent_[i0];
    785                     region2[iBack] = pivotValue * sign_[i0];
    786                     region[i0] = 0.0;
    787                     region[otherRow] += pivotValue;
    788                     i0 = otherRow;
    789                     double pivotValue1 = region[i1];
    790                     // put back now ?
    791                     int iBack1 = permuteBack_[i1];
    792                     regionIndex2[numberNonZero++] = iBack1;
    793                     int otherRow1 = parent_[i1];
    794                     region2[iBack1] = pivotValue1 * sign_[i1];
    795                     region[i1] = 0.0;
    796                     region[otherRow1] += pivotValue1;
    797                     i1 = otherRow1;
    798                }
    799           } else {
    800                // set up linked lists at each depth
    801                // stack2 is start, stack is next
    802                int greatestDepth = -1;
    803                //mark_[numberRows_]=1;
    804                for (i = 0; i < numberNonZero; i++) {
    805                     int j = regionIndex2[i];
    806                     double value = region2[j];
    807                     region2[j] = 0.0;
    808                     region[j] = value;
    809                     regionIndex[i] = j;
    810                     int iDepth = depth_[j];
    811                     if (iDepth > greatestDepth)
    812                          greatestDepth = iDepth;
    813                     // and back until marked
    814                     while (!mark_[j]) {
    815                          int iNext = stack2_[iDepth];
    816                          stack2_[iDepth] = j;
    817                          stack_[j] = iNext;
    818                          mark_[j] = 1;
    819                          iDepth--;
    820                          j = parent_[j];
    821                     }
    822                }
    823                numberNonZero = 0;
    824                for (; greatestDepth >= 0; greatestDepth--) {
    825                     int iPivot = stack2_[greatestDepth];
    826                     stack2_[greatestDepth] = -1;
    827                     while (iPivot >= 0) {
    828                          mark_[iPivot] = 0;
    829                          double pivotValue = region[iPivot];
    830                          if (pivotValue) {
    831                               // put back now ?
    832                               int iBack = permuteBack_[iPivot];
    833                               regionIndex2[numberNonZero++] = iBack;
    834                               int otherRow = parent_[iPivot];
    835                               region2[iBack] = pivotValue * sign_[iPivot];
    836                               region[iPivot] = 0.0;
    837                               region[otherRow] += pivotValue;
    838                          }
    839                          iPivot = stack_[iPivot];
    840                     }
    841                }
     721        }
     722      } else {
     723        for (; greatestDepth >= 0; greatestDepth--) {
     724          int iPivot = stack2_[greatestDepth];
     725          stack2_[greatestDepth] = -1;
     726          while (iPivot >= 0) {
     727            mark_[iPivot] = 0;
     728            double pivotValue = region[iPivot];
     729            if (pivotValue) {
     730              // put back now ?
     731              int iBack = permuteBack_[iPivot];
     732              double value = pivotValue * sign_[iPivot];
     733              region2[numberNonZero] = value;
     734              regionIndex2[numberNonZero++] = iBack;
     735              if (iBack == pivotRow)
     736                returnValue = value;
     737              int otherRow = parent_[iPivot];
     738              region[iPivot] = 0.0;
     739              region[otherRow] += pivotValue;
     740            }
     741            iPivot = stack_[iPivot];
    842742          }
    843           if (pivotRow >= 0)
    844                returnValue = region2[pivotRow];
    845      }
    846      region[numberRows_] = 0.0;
    847      regionSparse2->setNumElements(numberNonZero);
     743        }
     744      }
     745    }
     746  } else {
     747    if (doTwo && region2[i0] * region2[i1] < 0.0) {
     748      // If just +- 1 then could go backwards on depth until join
     749      region[i0] = region2[i0];
     750      region2[i0] = 0.0;
     751      region[i1] = region2[i1];
     752      region2[i1] = 0.0;
     753      int iDepth0 = depth_[i0];
     754      int iDepth1 = depth_[i1];
     755      if (iDepth1 > iDepth0) {
     756        int temp = i0;
     757        i0 = i1;
     758        i1 = temp;
     759        temp = iDepth0;
     760        iDepth0 = iDepth1;
     761        iDepth1 = temp;
     762      }
     763      numberNonZero = 0;
     764      while (iDepth0 > iDepth1) {
     765        double pivotValue = region[i0];
     766        // put back now ?
     767        int iBack = permuteBack_[i0];
     768        regionIndex2[numberNonZero++] = iBack;
     769        int otherRow = parent_[i0];
     770        region2[iBack] = pivotValue * sign_[i0];
     771        region[i0] = 0.0;
     772        region[otherRow] += pivotValue;
     773        iDepth0--;
     774        i0 = otherRow;
     775      }
     776      while (i0 != i1) {
     777        double pivotValue = region[i0];
     778        // put back now ?
     779        int iBack = permuteBack_[i0];
     780        regionIndex2[numberNonZero++] = iBack;
     781        int otherRow = parent_[i0];
     782        region2[iBack] = pivotValue * sign_[i0];
     783        region[i0] = 0.0;
     784        region[otherRow] += pivotValue;
     785        i0 = otherRow;
     786        double pivotValue1 = region[i1];
     787        // put back now ?
     788        int iBack1 = permuteBack_[i1];
     789        regionIndex2[numberNonZero++] = iBack1;
     790        int otherRow1 = parent_[i1];
     791        region2[iBack1] = pivotValue1 * sign_[i1];
     792        region[i1] = 0.0;
     793        region[otherRow1] += pivotValue1;
     794        i1 = otherRow1;
     795      }
     796    } else {
     797      // set up linked lists at each depth
     798      // stack2 is start, stack is next
     799      int greatestDepth = -1;
     800      //mark_[numberRows_]=1;
     801      for (i = 0; i < numberNonZero; i++) {
     802        int j = regionIndex2[i];
     803        double value = region2[j];
     804        region2[j] = 0.0;
     805        region[j] = value;
     806        regionIndex[i] = j;
     807        int iDepth = depth_[j];
     808        if (iDepth > greatestDepth)
     809          greatestDepth = iDepth;
     810        // and back until marked
     811        while (!mark_[j]) {
     812          int iNext = stack2_[iDepth];
     813          stack2_[iDepth] = j;
     814          stack_[j] = iNext;
     815          mark_[j] = 1;
     816          iDepth--;
     817          j = parent_[j];
     818        }
     819      }
     820      numberNonZero = 0;
     821      for (; greatestDepth >= 0; greatestDepth--) {
     822        int iPivot = stack2_[greatestDepth];
     823        stack2_[greatestDepth] = -1;
     824        while (iPivot >= 0) {
     825          mark_[iPivot] = 0;
     826          double pivotValue = region[iPivot];
     827          if (pivotValue) {
     828            // put back now ?
     829            int iBack = permuteBack_[iPivot];
     830            regionIndex2[numberNonZero++] = iBack;
     831            int otherRow = parent_[iPivot];
     832            region2[iBack] = pivotValue * sign_[iPivot];
     833            region[iPivot] = 0.0;
     834            region[otherRow] += pivotValue;
     835          }
     836          iPivot = stack_[iPivot];
     837        }
     838      }
     839    }
     840    if (pivotRow >= 0)
     841      returnValue = region2[pivotRow];
     842  }
     843  region[numberRows_] = 0.0;
     844  regionSparse2->setNumElements(numberNonZero);
    848845#ifdef FULL_DEBUG
    849      {
    850           int i;
    851           for (i = 0; i < numberRows_; i++) {
    852                assert(!mark_[i]);
    853                assert (stack2_[i] == -1);
    854           }
    855      }
     846  {
     847    int i;
     848    for (i = 0; i < numberRows_; i++) {
     849      assert(!mark_[i]);
     850      assert(stack2_[i] == -1);
     851    }
     852  }
    856853#endif
    857      return returnValue;
     854  return returnValue;
    858855}
    859856/* Updates one column (FTRAN) to/from array
     
    862859   have got code working using this simple method - thank you!
    863860   (the only exception is if you know input is dense e.g. rhs) */
    864 int
    865 ClpNetworkBasis::updateColumn (  CoinIndexedVector * regionSparse,
    866                                  double region2[] ) const
     861int ClpNetworkBasis::updateColumn(CoinIndexedVector *regionSparse,
     862  double region2[]) const
    867863{
    868      regionSparse->clear (  );
    869      double *region = regionSparse->denseVector (  );
    870      int numberNonZero = 0;
    871      int *regionIndex = regionSparse->getIndices (  );
    872      int i;
    873      // set up linked lists at each depth
    874      // stack2 is start, stack is next
    875      int greatestDepth = -1;
    876      for (i = 0; i < numberRows_; i++) {
    877           double value = region2[i];
    878           if (value) {
    879                region2[i] = 0.0;
    880                region[i] = value;
    881                regionIndex[numberNonZero++] = i;
    882                int j = i;
    883                int iDepth = depth_[j];
    884                if (iDepth > greatestDepth)
    885                     greatestDepth = iDepth;
    886                // and back until marked
    887                while (!mark_[j]) {
    888                     int iNext = stack2_[iDepth];
    889                     stack2_[iDepth] = j;
    890                     stack_[j] = iNext;
    891                     mark_[j] = 1;
    892                     iDepth--;
    893                     j = parent_[j];
    894                }
    895           }
    896      }
    897      numberNonZero = 0;
    898      for (; greatestDepth >= 0; greatestDepth--) {
    899           int iPivot = stack2_[greatestDepth];
    900           stack2_[greatestDepth] = -1;
    901           while (iPivot >= 0) {
    902                mark_[iPivot] = 0;
    903                double pivotValue = region[iPivot];
    904                if (pivotValue) {
    905                     // put back now ?
    906                     int iBack = permuteBack_[iPivot];
    907                     numberNonZero++;
    908                     int otherRow = parent_[iPivot];
    909                     region2[iBack] = pivotValue * sign_[iPivot];
    910                     region[iPivot] = 0.0;
    911                     region[otherRow] += pivotValue;
    912                }
    913                iPivot = stack_[iPivot];
    914           }
    915      }
    916      region[numberRows_] = 0.0;
    917      return numberNonZero;
     864  regionSparse->clear();
     865  double *region = regionSparse->denseVector();
     866  int numberNonZero = 0;
     867  int *regionIndex = regionSparse->getIndices();
     868  int i;
     869  // set up linked lists at each depth
     870  // stack2 is start, stack is next
     871  int greatestDepth = -1;
     872  for (i = 0; i < numberRows_; i++) {
     873    double value = region2[i];
     874    if (value) {
     875      region2[i] = 0.0;
     876      region[i] = value;
     877      regionIndex[numberNonZero++] = i;
     878      int j = i;
     879      int iDepth = depth_[j];
     880      if (iDepth > greatestDepth)
     881        greatestDepth = iDepth;
     882      // and back until marked
     883      while (!mark_[j]) {
     884        int iNext = stack2_[iDepth];
     885        stack2_[iDepth] = j;
     886        stack_[j] = iNext;
     887        mark_[j] = 1;
     888        iDepth--;
     889        j = parent_[j];
     890      }
     891    }
     892  }
     893  numberNonZero = 0;
     894  for (; greatestDepth >= 0; greatestDepth--) {
     895    int iPivot = stack2_[greatestDepth];
     896    stack2_[greatestDepth] = -1;
     897    while (iPivot >= 0) {
     898      mark_[iPivot] = 0;
     899      double pivotValue = region[iPivot];
     900      if (pivotValue) {
     901        // put back now ?
     902        int iBack = permuteBack_[iPivot];
     903        numberNonZero++;
     904        int otherRow = parent_[iPivot];
     905        region2[iBack] = pivotValue * sign_[iPivot];
     906        region[iPivot] = 0.0;
     907        region[otherRow] += pivotValue;
     908      }
     909      iPivot = stack_[iPivot];
     910    }
     911  }
     912  region[numberRows_] = 0.0;
     913  return numberNonZero;
    918914}
    919915/* Updates one column transpose (BTRAN)
     
    923919   (the only exception is if you know input is dense e.g. dense objective)
    924920   returns number of nonzeros */
    925 int
    926 ClpNetworkBasis::updateColumnTranspose (  CoinIndexedVector * regionSparse,
    927           double region2[] ) const
     921int ClpNetworkBasis::updateColumnTranspose(CoinIndexedVector *regionSparse,
     922  double region2[]) const
    928923{
    929      // permute in after copying
    930      // so will end up in right place
    931      double *region = regionSparse->denseVector (  );
    932      int *regionIndex = regionSparse->getIndices (  );
    933      int i;
    934      int numberNonZero = 0;
    935      CoinMemcpyN(region2, numberRows_, region);
    936      for (i = 0; i < numberRows_; i++) {
    937           double value = region[i];
    938           if (value) {
    939                int k = permute_[i];
    940                region[i] = 0.0;
    941                region2[k] = value;
    942                regionIndex[numberNonZero++] = k;
    943                mark_[k] = 1;
    944           }
    945      }
    946      // copy back
    947      // set up linked lists at each depth
    948      // stack2 is start, stack is next
    949      int greatestDepth = -1;
    950      int smallestDepth = numberRows_;
    951      for (i = 0; i < numberNonZero; i++) {
    952           int j = regionIndex[i];
    953           // add in
    954           int iDepth = depth_[j];
    955           smallestDepth = CoinMin(iDepth, smallestDepth) ;
    956           greatestDepth = CoinMax(iDepth, greatestDepth) ;
    957           int jNext = stack2_[iDepth];
    958           stack2_[iDepth] = j;
    959           stack_[j] = jNext;
    960           // and put all descendants on list
    961           int iChild = descendant_[j];
    962           while (iChild >= 0) {
    963                if (!mark_[iChild]) {
    964                     regionIndex[numberNonZero++] = iChild;
    965                     mark_[iChild] = 1;
    966                }
    967                iChild = rightSibling_[iChild];
    968           }
    969      }
    970      numberNonZero = 0;
    971      region2[numberRows_] = 0.0;
    972      int iDepth;
    973      for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
    974           int iPivot = stack2_[iDepth];
    975           stack2_[iDepth] = -1;
    976           while (iPivot >= 0) {
    977                mark_[iPivot] = 0;
    978                double pivotValue = region2[iPivot];
    979                int otherRow = parent_[iPivot];
    980                double otherValue = region2[otherRow];
    981                pivotValue = sign_[iPivot] * pivotValue + otherValue;
    982                region2[iPivot] = pivotValue;
    983                if (pivotValue)
    984                     numberNonZero++;
    985                iPivot = stack_[iPivot];
    986           }
    987      }
    988      return numberNonZero;
     924  // permute in after copying
     925  // so will end up in right place
     926  double *region = regionSparse->denseVector();
     927  int *regionIndex = regionSparse->getIndices();
     928  int i;
     929  int numberNonZero = 0;
     930  CoinMemcpyN(region2, numberRows_, region);
     931  for (i = 0; i < numberRows_; i++) {
     932    double value = region[i];
     933    if (value) {
     934      int k = permute_[i];
     935      region[i] = 0.0;
     936      region2[k] = value;
     937      regionIndex[numberNonZero++] = k;
     938      mark_[k] = 1;
     939    }
     940  }
     941  // copy back
     942  // set up linked lists at each depth
     943  // stack2 is start, stack is next
     944  int greatestDepth = -1;
     945  int smallestDepth = numberRows_;
     946  for (i = 0; i < numberNonZero; i++) {
     947    int j = regionIndex[i];
     948    // add in
     949    int iDepth = depth_[j];
     950    smallestDepth = CoinMin(iDepth, smallestDepth);
     951    greatestDepth = CoinMax(iDepth, greatestDepth);
     952    int jNext = stack2_[iDepth];
     953    stack2_[iDepth] = j;
     954    stack_[j] = jNext;
     955    // and put all descendants on list
     956    int iChild = descendant_[j];
     957    while (iChild >= 0) {
     958      if (!mark_[iChild]) {
     959        regionIndex[numberNonZero++] = iChild;
     960        mark_[iChild] = 1;
     961      }
     962      iChild = rightSibling_[iChild];
     963    }
     964  }
     965  numberNonZero = 0;
     966  region2[numberRows_] = 0.0;
     967  int iDepth;
     968  for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
     969    int iPivot = stack2_[iDepth];
     970    stack2_[iDepth] = -1;
     971    while (iPivot >= 0) {
     972      mark_[iPivot] = 0;
     973      double pivotValue = region2[iPivot];
     974      int otherRow = parent_[iPivot];
     975      double otherValue = region2[otherRow];
     976      pivotValue = sign_[iPivot] * pivotValue + otherValue;
     977      region2[iPivot] = pivotValue;
     978      if (pivotValue)
     979        numberNonZero++;
     980      iPivot = stack_[iPivot];
     981    }
     982  }
     983  return numberNonZero;
    989984}
    990985/* Updates one column (BTRAN) from region2 */
    991 int
    992 ClpNetworkBasis::updateColumnTranspose (  CoinIndexedVector * regionSparse,
    993           CoinIndexedVector * regionSparse2) const
     986int ClpNetworkBasis::updateColumnTranspose(CoinIndexedVector *regionSparse,
     987  CoinIndexedVector *regionSparse2) const
    994988{
    995      // permute in - presume small number so copy back
    996      // so will end up in right place
    997      regionSparse->clear (  );
    998      double *region = regionSparse->denseVector (  );
    999      double *region2 = regionSparse2->denseVector (  );
    1000      int *regionIndex2 = regionSparse2->getIndices (  );
    1001      int numberNonZero2 = regionSparse2->getNumElements (  );
    1002      int *regionIndex = regionSparse->getIndices (  );
    1003      int i;
    1004      int numberNonZero = 0;
    1005      bool packed = regionSparse2->packedMode();
    1006      if (packed) {
    1007           for (i = 0; i < numberNonZero2; i++) {
    1008                int k = regionIndex2[i];
    1009                int j = permute_[k];
    1010                double value = region2[i];
    1011                region2[i] = 0.0;
    1012                region[j] = value;
    1013                mark_[j] = 1;
    1014                regionIndex[numberNonZero++] = j;
    1015           }
    1016           // set up linked lists at each depth
    1017           // stack2 is start, stack is next
    1018           int greatestDepth = -1;
    1019           int smallestDepth = numberRows_;
    1020           //mark_[numberRows_]=1;
    1021           for (i = 0; i < numberNonZero2; i++) {
    1022                int j = regionIndex[i];
    1023                regionIndex2[i] = j;
    1024                // add in
    1025                int iDepth = depth_[j];
    1026                smallestDepth = CoinMin(iDepth, smallestDepth) ;
    1027                greatestDepth = CoinMax(iDepth, greatestDepth) ;
    1028                int jNext = stack2_[iDepth];
    1029                stack2_[iDepth] = j;
    1030                stack_[j] = jNext;
    1031                // and put all descendants on list
    1032                int iChild = descendant_[j];
    1033                while (iChild >= 0) {
    1034                     if (!mark_[iChild]) {
    1035                          regionIndex2[numberNonZero++] = iChild;
    1036                          mark_[iChild] = 1;
    1037                     }
    1038                     iChild = rightSibling_[iChild];
    1039                }
    1040           }
    1041           for (; i < numberNonZero; i++) {
    1042                int j = regionIndex2[i];
    1043                // add in
    1044                int iDepth = depth_[j];
    1045                smallestDepth = CoinMin(iDepth, smallestDepth) ;
    1046                greatestDepth = CoinMax(iDepth, greatestDepth) ;
    1047                int jNext = stack2_[iDepth];
    1048                stack2_[iDepth] = j;
    1049                stack_[j] = jNext;
    1050                // and put all descendants on list
    1051                int iChild = descendant_[j];
    1052                while (iChild >= 0) {
    1053                     if (!mark_[iChild]) {
    1054                          regionIndex2[numberNonZero++] = iChild;
    1055                          mark_[iChild] = 1;
    1056                     }
    1057                     iChild = rightSibling_[iChild];
    1058                }
    1059           }
    1060           numberNonZero2 = 0;
    1061           region[numberRows_] = 0.0;
    1062           int iDepth;
    1063           for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
    1064                int iPivot = stack2_[iDepth];
    1065                stack2_[iDepth] = -1;
    1066                while (iPivot >= 0) {
    1067                     mark_[iPivot] = 0;
    1068                     double pivotValue = region[iPivot];
    1069                     int otherRow = parent_[iPivot];
    1070                     double otherValue = region[otherRow];
    1071                     pivotValue = sign_[iPivot] * pivotValue + otherValue;
    1072                     region[iPivot] = pivotValue;
    1073                     if (pivotValue) {
    1074                          region2[numberNonZero2] = pivotValue;
    1075                          regionIndex2[numberNonZero2++] = iPivot;
    1076                     }
    1077                     iPivot = stack_[iPivot];
    1078                }
    1079           }
    1080           // zero out
    1081           for (i = 0; i < numberNonZero2; i++) {
    1082                int k = regionIndex2[i];
    1083                region[k] = 0.0;
    1084           }
    1085      } else {
    1086           for (i = 0; i < numberNonZero2; i++) {
    1087                int k = regionIndex2[i];
    1088                int j = permute_[k];
    1089                double value = region2[k];
    1090                region2[k] = 0.0;
    1091                region[j] = value;
    1092                mark_[j] = 1;
    1093                regionIndex[numberNonZero++] = j;
    1094           }
    1095           // copy back
    1096           // set up linked lists at each depth
    1097           // stack2 is start, stack is next
    1098           int greatestDepth = -1;
    1099           int smallestDepth = numberRows_;
    1100           //mark_[numberRows_]=1;
    1101           for (i = 0; i < numberNonZero2; i++) {
    1102                int j = regionIndex[i];
    1103                double value = region[j];
    1104                region[j] = 0.0;
    1105                region2[j] = value;
    1106                regionIndex2[i] = j;
    1107                // add in
    1108                int iDepth = depth_[j];
    1109                smallestDepth = CoinMin(iDepth, smallestDepth) ;
    1110                greatestDepth = CoinMax(iDepth, greatestDepth) ;
    1111                int jNext = stack2_[iDepth];
    1112                stack2_[iDepth] = j;
    1113                stack_[j] = jNext;
    1114                // and put all descendants on list
    1115                int iChild = descendant_[j];
    1116                while (iChild >= 0) {
    1117                     if (!mark_[iChild]) {
    1118                          regionIndex2[numberNonZero++] = iChild;
    1119                          mark_[iChild] = 1;
    1120                     }
    1121                     iChild = rightSibling_[iChild];
    1122                }
    1123           }
    1124           for (; i < numberNonZero; i++) {
    1125                int j = regionIndex2[i];
    1126                // add in
    1127                int iDepth = depth_[j];
    1128                smallestDepth = CoinMin(iDepth, smallestDepth) ;
    1129                greatestDepth = CoinMax(iDepth, greatestDepth) ;
    1130                int jNext = stack2_[iDepth];
    1131                stack2_[iDepth] = j;
    1132                stack_[j] = jNext;
    1133                // and put all descendants on list
    1134                int iChild = descendant_[j];
    1135                while (iChild >= 0) {
    1136                     if (!mark_[iChild]) {
    1137                          regionIndex2[numberNonZero++] = iChild;
    1138                          mark_[iChild] = 1;
    1139                     }
    1140                     iChild = rightSibling_[iChild];
    1141                }
    1142           }
    1143           numberNonZero2 = 0;
    1144           region2[numberRows_] = 0.0;
    1145           int iDepth;
    1146           for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
    1147                int iPivot = stack2_[iDepth];
    1148                stack2_[iDepth] = -1;
    1149                while (iPivot >= 0) {
    1150                     mark_[iPivot] = 0;
    1151                     double pivotValue = region2[iPivot];
    1152                     int otherRow = parent_[iPivot];
    1153                     double otherValue = region2[otherRow];
    1154                     pivotValue = sign_[iPivot] * pivotValue + otherValue;
    1155                     region2[iPivot] = pivotValue;
    1156                     if (pivotValue)
    1157                          regionIndex2[numberNonZero2++] = iPivot;
    1158                     iPivot = stack_[iPivot];
    1159                }
    1160           }
    1161      }
    1162      regionSparse2->setNumElements(numberNonZero2);
     989  // permute in - presume small number so copy back
     990  // so will end up in right place
     991  regionSparse->clear();
     992  double *region = regionSparse->denseVector();
     993  double *region2 = regionSparse2->denseVector();
     994  int *regionIndex2 = regionSparse2->getIndices();
     995  int numberNonZero2 = regionSparse2->getNumElements();
     996  int *regionIndex = regionSparse->getIndices();
     997  int i;
     998  int numberNonZero = 0;
     999  bool packed = regionSparse2->packedMode();
     1000  if (packed) {
     1001    for (i = 0; i < numberNonZero2; i++) {
     1002      int k = regionIndex2[i];
     1003      int j = permute_[k];
     1004      double value = region2[i];
     1005      region2[i] = 0.0;
     1006      region[j] = value;
     1007      mark_[j] = 1;
     1008      regionIndex[numberNonZero++] = j;
     1009    }
     1010    // set up linked lists at each depth
     1011    // stack2 is start, stack is next
     1012    int greatestDepth = -1;
     1013    int smallestDepth = numberRows_;
     1014    //mark_[numberRows_]=1;
     1015    for (i = 0; i < numberNonZero2; i++) {
     1016      int j = regionIndex[i];
     1017      regionIndex2[i] = j;
     1018      // add in
     1019      int iDepth = depth_[j];
     1020      smallestDepth = CoinMin(iDepth, smallestDepth);
     1021      greatestDepth = CoinMax(iDepth, greatestDepth);
     1022      int jNext = stack2_[iDepth];
     1023      stack2_[iDepth] = j;
     1024      stack_[j] = jNext;
     1025      // and put all descendants on list
     1026      int iChild = descendant_[j];
     1027      while (iChild >= 0) {
     1028        if (!mark_[iChild]) {
     1029          regionIndex2[numberNonZero++] = iChild;
     1030          mark_[iChild] = 1;
     1031        }
     1032        iChild = rightSibling_[iChild];
     1033      }
     1034    }
     1035    for (; i < numberNonZero; i++) {
     1036      int j = regionIndex2[i];
     1037      // add in
     1038      int iDepth = depth_[j];
     1039      smallestDepth = CoinMin(iDepth, smallestDepth);
     1040      greatestDepth = CoinMax(iDepth, greatestDepth);
     1041      int jNext = stack2_[iDepth];
     1042      stack2_[iDepth] = j;
     1043      stack_[j] = jNext;
     1044      // and put all descendants on list
     1045      int iChild = descendant_[j];
     1046      while (iChild >= 0) {
     1047        if (!mark_[iChild]) {
     1048          regionIndex2[numberNonZero++] = iChild;
     1049          mark_[iChild] = 1;
     1050        }
     1051        iChild = rightSibling_[iChild];
     1052      }
     1053    }
     1054    numberNonZero2 = 0;
     1055    region[numberRows_] = 0.0;
     1056    int iDepth;
     1057    for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
     1058      int iPivot = stack2_[iDepth];
     1059      stack2_[iDepth] = -1;
     1060      while (iPivot >= 0) {
     1061        mark_[iPivot] = 0;
     1062        double pivotValue = region[iPivot];
     1063        int otherRow = parent_[iPivot];
     1064        double otherValue = region[otherRow];
     1065        pivotValue = sign_[iPivot] * pivotValue + otherValue;
     1066        region[iPivot] = pivotValue;
     1067        if (pivotValue) {
     1068          region2[numberNonZero2] = pivotValue;
     1069          regionIndex2[numberNonZero2++] = iPivot;
     1070        }
     1071        iPivot = stack_[iPivot];
     1072      }
     1073    }
     1074    // zero out
     1075    for (i = 0; i < numberNonZero2; i++) {
     1076      int k = regionIndex2[i];
     1077      region[k] = 0.0;
     1078    }
     1079  } else {
     1080    for (i = 0; i < numberNonZero2; i++) {
     1081      int k = regionIndex2[i];
     1082      int j = permute_[k];
     1083      double value = region2[k];
     1084      region2[k] = 0.0;
     1085      region[j] = value;
     1086      mark_[j] = 1;
     1087      regionIndex[numberNonZero++] = j;
     1088    }
     1089    // copy back
     1090    // set up linked lists at each depth
     1091    // stack2 is start, stack is next
     1092    int greatestDepth = -1;
     1093    int smallestDepth = numberRows_;
     1094    //mark_[numberRows_]=1;
     1095    for (i = 0; i < numberNonZero2; i++) {
     1096      int j = regionIndex[i];
     1097      double value = region[j];
     1098      region[j] = 0.0;
     1099      region2[j] = value;
     1100      regionIndex2[i] = j;
     1101      // add in
     1102      int iDepth = depth_[j];
     1103      smallestDepth = CoinMin(iDepth, smallestDepth);
     1104      greatestDepth = CoinMax(iDepth, greatestDepth);
     1105      int jNext = stack2_[iDepth];
     1106      stack2_[iDepth] = j;
     1107      stack_[j] = jNext;
     1108      // and put all descendants on list
     1109      int iChild = descendant_[j];
     1110      while (iChild >= 0) {
     1111        if (!mark_[iChild]) {
     1112          regionIndex2[numberNonZero++] = iChild;
     1113          mark_[iChild] = 1;
     1114        }
     1115        iChild = rightSibling_[iChild];
     1116      }
     1117    }
     1118    for (; i < numberNonZero; i++) {
     1119      int j = regionIndex2[i];
     1120      // add in
     1121      int iDepth = depth_[j];
     1122      smallestDepth = CoinMin(iDepth, smallestDepth);
     1123      greatestDepth = CoinMax(iDepth, greatestDepth);
     1124      int jNext = stack2_[iDepth];
     1125      stack2_[iDepth] = j;
     1126      stack_[j] = jNext;
     1127      // and put all descendants on list
     1128      int iChild = descendant_[j];
     1129      while (iChild >= 0) {
     1130        if (!mark_[iChild]) {
     1131          regionIndex2[numberNonZero++] = iChild;
     1132          mark_[iChild] = 1;
     1133        }
     1134        iChild = rightSibling_[iChild];
     1135      }
     1136    }
     1137    numberNonZero2 = 0;
     1138    region2[numberRows_] = 0.0;
     1139    int iDepth;
     1140    for (iDepth = smallestDepth; iDepth <= greatestDepth; iDepth++) {
     1141      int iPivot = stack2_[iDepth];
     1142      stack2_[iDepth] = -1;
     1143      while (iPivot >= 0) {
     1144        mark_[iPivot] = 0;
     1145        double pivotValue = region2[iPivot];
     1146        int otherRow = parent_[iPivot];
     1147        double otherValue = region2[otherRow];
     1148        pivotValue = sign_[iPivot] * pivotValue + otherValue;
     1149        region2[iPivot] = pivotValue;
     1150        if (pivotValue)
     1151          regionIndex2[numberNonZero2++] = iPivot;
     1152        iPivot = stack_[iPivot];
     1153      }
     1154    }
     1155  }
     1156  regionSparse2->setNumElements(numberNonZero2);
    11631157#ifdef FULL_DEBUG
    1164      {
    1165           int i;
    1166           for (i = 0; i < numberRows_; i++) {
    1167                assert(!mark_[i]);
    1168                assert (stack2_[i] == -1);
    1169           }
    1170      }
     1158  {
     1159    int i;
     1160    for (i = 0; i < numberRows_; i++) {
     1161      assert(!mark_[i]);
     1162      assert(stack2_[i] == -1);
     1163    }
     1164  }
    11711165#endif
    1172      return numberNonZero2;
     1166  return numberNonZero2;
    11731167}
     1168
     1169/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
     1170*/
Note: See TracChangeset for help on using the changeset viewer.