ldbl_ratio.c

Go to the documentation of this file.
00001 /****************************************************************************/
00002 /*                                                                          */
00003 /*  This file is part of QSopt_ex.                                          */
00004 /*                                                                          */
00005 /*  (c) Copyright 2006 by David Applegate, William Cook, Sanjeeb Dash,      */
00006 /*  and Daniel Espinoza                                                     */
00007 /*                                                                          */
00008 /*  Sanjeeb Dash ownership of copyright in QSopt_ex is derived from his     */
00009 /*  copyright in QSopt.                                                     */
00010 /*                                                                          */
00011 /*  This code may be used under the terms of the GNU General Public License */
00012 /*  (Version 2.1 or later) as published by the Free Software Foundation.    */
00013 /*                                                                          */
00014 /*  Alternatively, use is granted for research purposes only.               */
00015 /*                                                                          */
00016 /*  It is your choice of which of these two licenses you are operating      */
00017 /*  under.                                                                  */
00018 /*                                                                          */
00019 /*  We make no guarantees about the correctness or usefulness of this code. */
00020 /*                                                                          */
00021 /****************************************************************************/
00022 
00023 /* "$RCSfile: ldbl_ratio.c,v $ $Revision: 1.2 $ $Date: 2003/11/05 16:49:52 $"; */
00024 static int TRACE = 0;
00025 
00026 #include "config.h"
00027 #include "ldbl_sortrus.h"
00028 #include "stddefs.h"
00029 #include "ldbl_iqsutil.h"
00030 #include "ldbl_lpdefs.h"
00031 #include "ldbl_ratio.h"
00032 #include "ldbl_fct.h"
00033 #ifdef USEDMALLOC
00034 #include "dmalloc.h"
00035 #endif
00036 
00037 void ldbl_ILLratio_pI_test (
00038   ldbl_lpinfo * lp,
00039   int eindex,
00040   int dir,
00041   ldbl_ratio_res * rs)
00042 {
00043   int i = 0, k = 0;
00044   int col, ecol;
00045   int cbnd, indx = 0;
00046   int tctr = 0;
00047   int *perm = lp->upd.perm;
00048   int *ix = lp->upd.ix;
00049   long double *pivtol = &(lp->tol->pivot_tol);
00050   long double *dftol = &(lp->tol->id_tol);
00051 
00052    /*HHH*/ long double * t = lp->upd.t;
00053   long double t_i, delta, y_ij, rcost, nrcost, ntmp;
00054   long double *x, *l, *u;
00055 
00056    /*HHH*/ ldbl_EGlpNumInitVar (t_i);
00057   ldbl_EGlpNumInitVar (delta);
00058   ldbl_EGlpNumInitVar (y_ij);
00059   ldbl_EGlpNumInitVar (rcost);
00060   ldbl_EGlpNumInitVar (nrcost);
00061   ldbl_EGlpNumInitVar (ntmp);
00062   ldbl_EGlpNumZero (t_i);
00063   ldbl_EGlpNumZero (y_ij);
00064   ldbl_EGlpNumZero (delta);
00065   rs->lindex = -1;
00066   ldbl_EGlpNumZero (rs->tz);
00067   ldbl_EGlpNumZero (rs->pivotval);
00068   rs->ratio_stat = RATIO_FAILED;
00069   rs->lvstat = -1;
00070   ecol = lp->nbaz[eindex];
00071   ILL_IFTRACE2 ("%s:%d:%d:%d:%d", __func__, eindex, dir, ecol,
00072                 (VBOUNDED == lp->vtype[ecol]));
00073   if (lp->vtype[ecol] == VBOUNDED)
00074   {
00075     ldbl_EGlpNumCopyDiff (t[0], lp->uz[ecol], lp->lz[ecol]);
00076     ix[0] = BBOUND;
00077     ILL_IFTRACE2 (":%d[%d](%la,%la,%la)\n", ix[tctr], tctr,
00078                   ldbl_EGlpNumToLf (t[tctr]), ldbl_EGlpNumToLf (lp->uz[ecol]),
00079                   ldbl_EGlpNumToLf (lp->lz[ecol]));
00080     tctr++;
00081   }
00082   ILL_IFTRACE2 (":%d", lp->yjz.nzcnt);
00083   for (k = 0; k < lp->yjz.nzcnt; k++)
00084   {
00085     ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
00086     if (!ldbl_EGlpNumIsNeqZero (y_ij, *pivtol))
00087       continue;
00088 
00089     i = lp->yjz.indx[k];
00090     x = &(lp->xbz[i]);
00091     col = lp->baz[i];
00092     l = &(lp->lz[col]);
00093     u = &(lp->uz[col]);
00094 
00095     if ((dir == VINCREASE && ldbl_EGlpNumIsGreatZero (y_ij)) ||
00096         (dir == VDECREASE && ldbl_EGlpNumIsLessZero (y_ij)))
00097     {
00098       if (ldbl_EGlpNumIsLessZero (y_ij))
00099         ldbl_EGlpNumSign (y_ij);
00100       ILL_IFTRACE2 (":%d", lp->bfeas[i]);
00101       if (lp->bfeas[i] > 0)
00102       {
00103         ldbl_EGlpNumCopyDiffRatio (t[tctr], *x, *u, y_ij);
00104         ix[tctr] = 10 * k + BATOUPPER;
00105         ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr, ldbl_EGlpNumToLf (t[tctr]));
00106         tctr++;
00107         if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
00108         {
00109           ldbl_EGlpNumCopyDiffRatio (t[tctr], *x, *l, y_ij);
00110           ix[tctr] = 10 * k + BATOLOWER;
00111           ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr,
00112                         ldbl_EGlpNumToLf (t[tctr]));
00113           tctr++;
00114         }
00115       }
00116       else if (lp->bfeas[i] == 0)
00117       {
00118         if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
00119         {
00120           ldbl_EGlpNumCopyDiffRatio (t[tctr], *x, *l, y_ij);
00121           ix[tctr] = 10 * k + BATOLOWER;
00122           ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr,
00123                         ldbl_EGlpNumToLf (t[tctr]));
00124           tctr++;
00125         }
00126       }
00127     }
00128     else if ((dir == VINCREASE && ldbl_EGlpNumIsLessZero (y_ij)) ||
00129              (dir == VDECREASE && ldbl_EGlpNumIsGreatZero (y_ij)))
00130     {
00131       if (ldbl_EGlpNumIsLessZero (y_ij))
00132         ldbl_EGlpNumSign (y_ij);
00133       ILL_IFTRACE2 (":%d", lp->bfeas[i]);
00134       if (lp->bfeas[i] < 0)
00135       {
00136         ldbl_EGlpNumCopyDiffRatio (t[tctr], *l, *x, y_ij);
00137         ix[tctr] = 10 * k + BBTOLOWER;
00138         ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr, ldbl_EGlpNumToLf (t[tctr]));
00139         tctr++;
00140         if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
00141         {
00142           ldbl_EGlpNumCopyDiffRatio (t[tctr], *u, *x, y_ij);
00143           ix[tctr] = 10 * k + BBTOUPPER;
00144           ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr,
00145                         ldbl_EGlpNumToLf (t[tctr]));
00146           tctr++;
00147         }
00148       }
00149       else if (lp->bfeas[i] == 0)
00150       {
00151         if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
00152         {
00153           ldbl_EGlpNumCopyDiffRatio (t[tctr], *u, *x, y_ij);
00154           ix[tctr] = 10 * k + BBTOUPPER;
00155           ILL_IFTRACE2 (":%d[%d](%la)\n", ix[tctr], tctr,
00156                         ldbl_EGlpNumToLf (t[tctr]));
00157           tctr++;
00158         }
00159       }
00160     }
00161   }
00162   if (tctr == 0)
00163   {
00164     rs->ratio_stat = RATIO_FAILED;
00165     ILL_CLEANUP;
00166   }
00167 
00168   for (i = 0; i < tctr; i++)
00169     perm[i] = i;
00170   ldbl_ILLutil_EGlpNum_perm_quicksort (perm, t, tctr);
00171 
00172   ldbl_EGlpNumZero (lp->upd.c_obj);
00173   ldbl_EGlpNumCopy (rcost, lp->pIdz[eindex]);
00174   ILL_IFTRACE2 ("\n%s:%d:%lf", __func__, tctr, ldbl_EGlpNumToLf (rcost));
00175   for (i = 0; i < tctr; i++)
00176   {
00177     ldbl_EGlpNumCopy (t_i, t[perm[i]]);
00178     ldbl_EGlpNumCopy (ntmp, t_i);
00179     ldbl_EGlpNumSubTo (ntmp, delta);
00180     ldbl_EGlpNumAddInnProdTo (lp->upd.c_obj, ntmp, rcost);
00181     ldbl_EGlpNumCopy (delta, t_i);
00182     ILL_IFTRACE2 (":%d:%lf", perm[i], ldbl_EGlpNumToLf (delta));
00183      /*HHH*/ cbnd = ix[perm[i]] % 10;
00184     if (cbnd != BBOUND)
00185     {
00186       k = ix[perm[i]] / 10;
00187       ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
00188       indx = lp->yjz.indx[k];
00189       ILL_IFTRACE2 (":%d", indx);
00190     }
00191 
00192     switch (cbnd)
00193     {
00194     case BBOUND:
00195       rs->ratio_stat = RATIO_NOBCHANGE;
00196       ldbl_EGlpNumCopy (rs->tz, t_i);
00197       if (dir != VINCREASE)
00198         ldbl_EGlpNumSign (rs->tz);
00199       ILL_CLEANUP;
00200 
00201     case BATOLOWER:
00202     case BATOUPPER:
00203       ldbl_EGlpNumAddTo (rcost, y_ij);
00204       break;
00205     case BBTOLOWER:
00206     case BBTOUPPER:
00207       ldbl_EGlpNumSubTo (rcost, y_ij);
00208       break;
00209     }
00210     ldbl_EGlpNumCopyNeg (nrcost, rcost);
00211     if ((dir == VINCREASE && ldbl_EGlpNumIsLeq (nrcost, *dftol)) ||
00212         (dir == VDECREASE && ldbl_EGlpNumIsLeq (rcost, *dftol)))
00213     {
00214       /* change 5 to -1 if t_i > 0 is required below */
00215       if (ldbl_EGlpNumIsLessZero (t_i) && i > 5)
00216       {
00217         /* printf ("pIhell %.5f %d\n", t_i, i); */
00218         ldbl_EGlpNumDivUiTo (t_i, 2);
00219         rs->ratio_stat = RATIO_NEGATIVE;
00220         ldbl_EGlpNumZero (rs->tz);
00221         ILL_CLEANUP;
00222       }
00223       rs->lindex = indx;
00224       rs->ratio_stat = RATIO_BCHANGE;
00225       if (cbnd == BATOLOWER || cbnd == BBTOLOWER)
00226         rs->lvstat = STAT_LOWER;
00227       else
00228         rs->lvstat = STAT_UPPER;
00229 
00230       ldbl_EGlpNumCopy (rs->pivotval, y_ij);
00231       ldbl_EGlpNumCopy (rs->tz, t_i);
00232       if (dir != VINCREASE)
00233         ldbl_EGlpNumSign (rs->tz);
00234       ILL_CLEANUP;
00235     }
00236   }
00237 
00238 CLEANUP:
00239   ldbl_ILLfct_update_counts (lp, CNT_PIPIV, 0, rs->pivotval);
00240   ILL_IFTRACE2 (":tctr %d:%d\n", tctr, rs->ratio_stat);
00241   lp->upd.tctr = tctr;
00242   lp->upd.i = i;
00243   ldbl_EGlpNumCopy (lp->upd.tz, t_i);
00244   ldbl_EGlpNumCopy (lp->upd.piv, rs->pivotval);
00245   if (dir == VDECREASE)
00246     ldbl_EGlpNumSign (lp->upd.c_obj);
00247   if (rs->lindex != -1)
00248     lp->upd.fs = lp->bfeas[rs->lindex];
00249   ldbl_EGlpNumClearVar (t_i);
00250   ldbl_EGlpNumClearVar (delta);
00251   ldbl_EGlpNumClearVar (y_ij);
00252   ldbl_EGlpNumClearVar (rcost);
00253   ldbl_EGlpNumClearVar (nrcost);
00254   ldbl_EGlpNumClearVar (ntmp);
00255 }
00256 
00257 void ldbl_ILLratio_pII_test (
00258   ldbl_lpinfo * lp,
00259   int eindex,
00260   int dir,
00261   ldbl_ratio_res * rs)
00262 {
00263   int i, k, indx, col, ecol;
00264   long double *x, *l, *u, t_max, ayi_max, yi_max, ay_ij, y_ij, t_i, t_z;
00265   long double *pivtol = &(lp->tol->pivot_tol);
00266   long double *pftol = &(lp->tol->pfeas_tol);
00267 
00268   ldbl_EGlpNumInitVar (y_ij);
00269   ldbl_EGlpNumInitVar (ay_ij);
00270   ldbl_EGlpNumInitVar (t_i);
00271   ldbl_EGlpNumInitVar (t_z);
00272   ldbl_EGlpNumInitVar (t_max);
00273   ldbl_EGlpNumInitVar (yi_max);
00274   ldbl_EGlpNumInitVar (ayi_max);
00275    /*HHH*/ rs->boundch = 0;
00276   rs->lindex = -1;
00277   ldbl_EGlpNumZero (rs->tz);
00278   rs->ratio_stat = RATIO_FAILED;
00279   rs->lvstat = -1;
00280   ldbl_EGlpNumZero (rs->pivotval);
00281   ldbl_EGlpNumZero (rs->lbound);
00282   ecol = lp->nbaz[eindex];
00283 
00284   for (k = 0, ldbl_EGlpNumCopy (t_max, ldbl_INFTY); k < lp->yjz.nzcnt; k++)
00285   {
00286     ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
00287     ldbl_EGlpNumCopyAbs (ay_ij, y_ij);
00288     if (!ldbl_EGlpNumIsNeqZero (y_ij, *pivtol))
00289       continue;
00290 
00291     ldbl_EGlpNumCopy (t_i, ldbl_INFTY);
00292     i = lp->yjz.indx[k];
00293     x = &(lp->xbz[i]);
00294     col = lp->baz[i];
00295     l = &(lp->lz[col]);
00296     u = &(lp->uz[col]);
00297 
00298     if ((dir == VINCREASE && ldbl_EGlpNumIsGreatZero (y_ij)) ||
00299         (dir == VDECREASE && ldbl_EGlpNumIsLessZero (y_ij)))
00300     {
00301       if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
00302       {
00303         ldbl_EGlpNumCopyDiff (t_i, *x, *l);
00304         ldbl_EGlpNumAddTo (t_i, *pftol);
00305         ldbl_EGlpNumDivTo (t_i, ay_ij);
00306       }
00307     }
00308     else if ((dir == VINCREASE && ldbl_EGlpNumIsLessZero (y_ij)) ||
00309              (dir == VDECREASE && ldbl_EGlpNumIsGreatZero (y_ij)))
00310     {
00311       if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
00312       {
00313         ldbl_EGlpNumCopySum (t_i, *u, *pftol);
00314         ldbl_EGlpNumSubTo (t_i, *x);
00315         ldbl_EGlpNumDivTo (t_i, ay_ij);
00316       }
00317     }
00318     if (ldbl_EGlpNumIsEqqual (t_i, ldbl_INFTY))
00319       continue;
00320 
00321     if (ldbl_EGlpNumIsLess (t_i, t_max))
00322     {
00323       /*HHH tind = i; yval = fabs (y_ij); tval = t_i - pftol/fabs(y_ij); */
00324       ldbl_EGlpNumCopy (t_max, t_i);
00325     }
00326   }
00327   /* we use yi_max as temporal variable here */
00328   ldbl_EGlpNumCopyDiff (yi_max, lp->uz[ecol], lp->lz[ecol]);
00329   if (lp->vtype[ecol] == VBOUNDED && ldbl_EGlpNumIsLeq (yi_max, t_max))
00330   {
00331 
00332     ldbl_EGlpNumCopy (t_max, yi_max);
00333     rs->ratio_stat = RATIO_NOBCHANGE;
00334     ldbl_EGlpNumCopy (rs->tz, t_max);
00335     if (dir != VINCREASE)
00336       ldbl_EGlpNumSign (rs->tz);
00337     ILL_CLEANUP;
00338   }
00339 
00340   if (ldbl_EGlpNumIsLeq (ldbl_INFTY, t_max))
00341   {
00342     rs->ratio_stat = RATIO_UNBOUNDED;
00343     ILL_CLEANUP;
00344   }
00345   /*if (ldbl_EGlpNumIsLess (t_max, ldbl_zeroLpNum))
00346    * printf ("pIIhell\n");
00347    */
00348   indx = -1;
00349   ldbl_EGlpNumZero (t_z);
00350   ldbl_EGlpNumZero (yi_max);
00351   ldbl_EGlpNumZero (ayi_max);
00352   ILL_IFTRACE2 (":%d", lp->yjz.nzcnt);
00353   for (k = 0; k < lp->yjz.nzcnt; k++)
00354   {
00355     ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
00356     ldbl_EGlpNumCopyAbs (ay_ij, y_ij);
00357     if (!ldbl_EGlpNumIsNeqZero (y_ij, *pivtol))
00358       continue;
00359 
00360     ldbl_EGlpNumCopy (t_i, ldbl_INFTY);
00361     i = lp->yjz.indx[k];
00362     x = &(lp->xbz[i]);
00363     col = lp->baz[i];
00364     l = &(lp->lz[col]);
00365     u = &(lp->uz[col]);
00366 
00367     if ((dir == VINCREASE && ldbl_EGlpNumIsGreatZero (y_ij)) ||
00368         (dir == VDECREASE && ldbl_EGlpNumIsLessZero (y_ij)))
00369     {
00370       if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
00371         ldbl_EGlpNumCopyDiffRatio (t_i, *x, *l, ay_ij);
00372     }
00373     else if ((dir == VINCREASE && ldbl_EGlpNumIsLessZero (y_ij)) ||
00374              (dir == VDECREASE && ldbl_EGlpNumIsGreatZero (y_ij)))
00375     {
00376       if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
00377         ldbl_EGlpNumCopyDiffRatio (t_i, *u, *x, ay_ij);
00378     }
00379 
00380     if (ldbl_EGlpNumIsLeq (t_i, t_max))
00381     {
00382       if (ldbl_EGlpNumIsLess (ayi_max, ay_ij))
00383       {
00384         ldbl_EGlpNumCopy (yi_max, y_ij);
00385         ldbl_EGlpNumCopy (ayi_max, ay_ij);
00386         indx = i;
00387         ldbl_EGlpNumCopy (t_z, t_i);
00388         ILL_IFTRACE2 (":%d:%lf:%lf:%lf:%lf", indx, ldbl_EGlpNumToLf (t_i),
00389                       ldbl_EGlpNumToLf (t_max), ldbl_EGlpNumToLf (ayi_max),
00390                       ldbl_EGlpNumToLf (ay_ij));
00391       }
00392     }
00393   }
00394 
00395   if (indx < 0)
00396   {
00397     rs->ratio_stat = RATIO_FAILED;
00398   }
00399   else
00400   {
00401     /*
00402      * if (tind != rs->lindex){
00403      * HHHprintf ("tmax %e tval = %e yval = %e tind = %d\n", t_max, tval, yval, tind);
00404      * HHHprintf ("h tval = %e yval = %e tind = %d\n",rs->tz, yi_max, rs->lindex);
00405      * }
00406      */
00407     ILL_IFTRACE2 (":%d", indx);
00408     rs->lindex = indx;
00409     ldbl_EGlpNumCopy (rs->tz, t_z);
00410     ldbl_EGlpNumCopy (rs->pivotval, yi_max);
00411     rs->ratio_stat = RATIO_BCHANGE;
00412 
00413     if (dir == VINCREASE)
00414       rs->lvstat =
00415         (ldbl_EGlpNumIsGreatZero (yi_max)) ? STAT_LOWER : STAT_UPPER;
00416     else
00417       rs->lvstat =
00418         (ldbl_EGlpNumIsGreatZero (yi_max)) ? STAT_UPPER : STAT_LOWER;
00419 
00420     if (ldbl_EGlpNumIsLessZero (rs->tz))
00421     {
00422       ILL_IFTRACE2 ("need to change bound, tz=%la\n", ldbl_EGlpNumToLf (rs->tz));
00423       ldbl_EGlpNumCopyAbs (rs->tz, t_max);
00424       ldbl_EGlpNumDivUiTo (rs->tz, 10);
00425       rs->boundch = 1;
00426       ldbl_EGlpNumCopy (rs->lbound, lp->xbz[rs->lindex]);
00427       if (rs->lvstat == STAT_LOWER)
00428         ldbl_EGlpNumSubInnProdTo (rs->lbound, rs->tz, ayi_max);
00429       else
00430         ldbl_EGlpNumAddInnProdTo (rs->lbound, rs->tz, ayi_max);
00431     }
00432     if (dir == VDECREASE)
00433       ldbl_EGlpNumSign (rs->tz);
00434   }
00435 CLEANUP:
00436   ldbl_ILLfct_update_counts (lp, CNT_PIIPIV, 0, rs->pivotval);
00437   ldbl_EGlpNumClearVar (y_ij);
00438   ldbl_EGlpNumClearVar (ay_ij);
00439   ldbl_EGlpNumClearVar (t_i);
00440   ldbl_EGlpNumClearVar (t_z);
00441   ldbl_EGlpNumClearVar (t_max);
00442   ldbl_EGlpNumClearVar (yi_max);
00443   ldbl_EGlpNumClearVar (ayi_max);
00444 }
00445 
00446 #define ldbl_GET_XY_DRATIOTEST \
00447       if (lp->vstat[col] == STAT_UPPER){ \
00448         ldbl_EGlpNumCopyNeg(x,lp->dz[j]);\
00449         ldbl_EGlpNumCopy(y, *zAj);\
00450       } \
00451       else{ \
00452          ldbl_EGlpNumCopy(x, lp->dz[j]); \
00453          ldbl_EGlpNumCopyNeg(y, *zAj);\
00454       } \
00455       if (lvstat == STAT_UPPER) \
00456          ldbl_EGlpNumSign(y);
00457 
00458 
00459 void ldbl_ILLratio_dI_test (
00460   ldbl_lpinfo * lp,
00461   int lindex,
00462   int lvstat,
00463   ldbl_ratio_res * rs)
00464 {
00465   int j = 0, k;
00466   int col;
00467   int cbnd, indx;
00468   int tctr = 0;
00469   int *perm = lp->upd.perm;
00470   int *ix = lp->upd.ix;
00471   long double *t = lp->upd.t;
00472   long double *zAj, x, y, t_j, theta, rcost, delta;
00473   long double *pftol = &(lp->tol->ip_tol);
00474   long double *pivtol = &(lp->tol->pivot_tol);
00475 
00476   ldbl_EGlpNumInitVar (x);
00477   ldbl_EGlpNumInitVar (y);
00478   ldbl_EGlpNumInitVar (t_j);
00479   ldbl_EGlpNumInitVar (theta);
00480   ldbl_EGlpNumInitVar (rcost);
00481   ldbl_EGlpNumInitVar (delta);
00482   ldbl_EGlpNumZero (delta);
00483   ldbl_EGlpNumZero (t_j);
00484   ldbl_EGlpNumZero (rs->tz);
00485    /*HHH*/ rs->eindex = -1;
00486   rs->ratio_stat = RATIO_FAILED;
00487   ldbl_EGlpNumZero (rs->pivotval);
00488 
00489   for (k = 0; k < lp->zA.nzcnt; k++)
00490   {
00491     zAj = &(lp->zA.coef[k]);
00492     if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
00493       continue;
00494 
00495     ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
00496     j = lp->zA.indx[k];
00497     col = lp->nbaz[j];
00498 
00499     if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED)
00500       continue;
00501 
00502     ldbl_GET_XY_DRATIOTEST;
00503 
00504     if (ldbl_EGlpNumIsLessZero (y))
00505     {
00506       if (lp->dfeas[j] != 0 && lp->vstat[col] != STAT_ZERO)
00507       {
00508         ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00509         ix[tctr] = 10 * k + BBTOLOWER;
00510         tctr++;
00511       }
00512       else if (lp->vstat[col] == STAT_ZERO)
00513       {
00514         if (lp->dfeas[j] < 0)
00515         {
00516           ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00517           ix[tctr] = 10 * k + BBTOLOWER;
00518           tctr++;
00519         }
00520         if (lp->dfeas[j] <= 0)
00521         {
00522           ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00523           ix[tctr] = 10 * k + BBTOUPPER;
00524           tctr++;
00525         }
00526       }
00527     }
00528     else
00529     {
00530       if (lp->dfeas[j] > 0)
00531       {
00532         if (lp->vstat[col] == STAT_ZERO)
00533         {
00534           ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00535           ix[tctr] = 10 * k + BATOUPPER;
00536           tctr++;
00537           ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00538           ix[tctr] = 10 * k + BATOLOWER;
00539           tctr++;
00540         }
00541       }
00542       else if (lp->dfeas[j] == 0)
00543       {
00544         ldbl_EGlpNumCopyFrac (t[tctr], x, y);
00545         if (lp->vtype[col] == VBOUNDED)
00546           ix[tctr] = 10 * k + BSKIP;
00547         else
00548           ix[tctr] = 10 * k + BATOLOWER;
00549         tctr++;
00550       }
00551     }
00552   }
00553 
00554   if (tctr == 0)
00555   {
00556     rs->ratio_stat = RATIO_FAILED;
00557     ILL_CLEANUP;
00558   }
00559 
00560   for (j = 0; j < tctr; j++)
00561     perm[j] = j;
00562   ldbl_ILLutil_EGlpNum_perm_quicksort (perm, t, tctr);
00563 
00564   ldbl_EGlpNumZero (lp->upd.c_obj);
00565   ldbl_EGlpNumCopy (rcost, lp->xbz[lindex]);
00566   if (lvstat == STAT_LOWER)
00567     ldbl_EGlpNumSign (rcost);
00568   for (j = 0; j < tctr; j++)
00569   {
00570     cbnd = ix[perm[j]] % 10;
00571     if (cbnd == BSKIP)
00572       continue;
00573 
00574     ldbl_EGlpNumCopy (t_j, t[perm[j]]);
00575     ldbl_EGlpNumCopy (x, t_j);
00576     ldbl_EGlpNumSubTo (x, delta);
00577     ldbl_EGlpNumAddInnProdTo (lp->upd.c_obj, x, rcost);
00578     ldbl_EGlpNumCopy (delta, t_j);
00579     k = ix[perm[j]] / 10;
00580     zAj = &(lp->zA.coef[k]);
00581     indx = lp->zA.indx[k];
00582 
00583     if (lp->vstat[lp->nbaz[indx]] == STAT_LOWER
00584         || lp->vstat[lp->nbaz[indx]] == STAT_ZERO)
00585       ldbl_EGlpNumCopyNeg (theta, *zAj);
00586     else
00587       ldbl_EGlpNumCopy (theta, *zAj);
00588 
00589     if (lvstat == STAT_UPPER)
00590       ldbl_EGlpNumSign (theta);
00591 
00592     switch (cbnd)
00593     {
00594     case BATOLOWER:
00595     case BATOUPPER:
00596       ldbl_EGlpNumSubTo (rcost, theta);
00597       break;
00598     case BBTOLOWER:
00599     case BBTOUPPER:
00600       ldbl_EGlpNumAddTo (rcost, theta);
00601       break;
00602     }
00603     if (ldbl_EGlpNumIsLeq (rcost, *pftol))
00604     {
00605       /* if (t_j < 0.0) printf ("dIhell\n"); */
00606       rs->eindex = indx;
00607       ldbl_EGlpNumCopy (rs->tz, t_j);
00608       ldbl_EGlpNumCopy (rs->pivotval, *zAj);
00609       rs->ratio_stat = RATIO_BCHANGE;
00610       ILL_CLEANUP;
00611     }
00612   }
00613 
00614 CLEANUP:
00615   ldbl_ILLfct_update_counts (lp, CNT_DIPIV, 0, rs->pivotval);
00616   ILL_IFTRACE2 ("%s:tctr %d\n", __func__, tctr);
00617   lp->upd.tctr = tctr;
00618   lp->upd.i = j;
00619   ldbl_EGlpNumCopyAbs (lp->upd.tz, t_j);
00620   ldbl_EGlpNumCopy (lp->upd.piv, rs->pivotval);
00621   if (rs->eindex != -1)
00622     lp->upd.fs = lp->dfeas[rs->eindex];
00623   ldbl_EGlpNumClearVar (x);
00624   ldbl_EGlpNumClearVar (y);
00625   ldbl_EGlpNumClearVar (t_j);
00626   ldbl_EGlpNumClearVar (theta);
00627   ldbl_EGlpNumClearVar (rcost);
00628   ldbl_EGlpNumClearVar (delta);
00629 }
00630 
00631 void ldbl_ILLratio_dII_test (
00632   ldbl_lpinfo * lp,
00633   /*int lindex,*/
00634   int lvstat,
00635   ldbl_ratio_res * rs)
00636 {
00637   int j, k, indx;
00638   int col, ecol;
00639   long double *zAj, azAj, az_max, x, y, t_j, z_max, t_max, t_z;
00640   long double *dftol = &(lp->tol->dfeas_tol);
00641   long double *pivtol = &(lp->tol->pivot_tol);
00642 
00643   ldbl_EGlpNumInitVar (x);
00644   ldbl_EGlpNumInitVar (y);
00645   ldbl_EGlpNumInitVar (t_j);
00646   ldbl_EGlpNumInitVar (z_max);
00647   ldbl_EGlpNumInitVar (t_max);
00648   ldbl_EGlpNumInitVar (az_max);
00649   ldbl_EGlpNumInitVar (azAj);
00650   ldbl_EGlpNumInitVar (t_z);
00651   ldbl_EGlpNumZero (t_j);
00652   rs->coeffch = 0;
00653   ldbl_EGlpNumZero (rs->ecoeff);
00654   rs->eindex = -1;
00655   rs->ratio_stat = RATIO_FAILED;
00656   ILL_IFTRACE2 ("%s:tctr %d\n", __func__, 0);
00657   lp->upd.tctr = 0;
00658   ldbl_EGlpNumZero (lp->upd.dty);
00659   for (k = 0, ldbl_EGlpNumCopy (t_max, ldbl_INFTY); k < lp->zA.nzcnt; k++)
00660   {
00661     zAj = &(lp->zA.coef[k]);
00662     if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
00663       continue;
00664 
00665     ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
00666     j = lp->zA.indx[k];
00667     col = lp->nbaz[j];
00668 
00669     if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED)
00670       continue;
00671 
00672     ldbl_GET_XY_DRATIOTEST;
00673 
00674 //#warning adding/substracting tolerances to used value, is it rigght?
00675     if (ldbl_EGlpNumIsGreatZero (y))
00676     {
00677       //t_j = (x + dftol) / y;
00678       ldbl_EGlpNumCopySum (t_j, x, *dftol);
00679       ldbl_EGlpNumDivTo (t_j, y);
00680     }
00681     else
00682     {
00683 //#warning adding/substracting tolerances to used value, is it rigght?
00684       if (lp->vstat[col] == STAT_ZERO)
00685         ldbl_EGlpNumCopyDiffRatio (t_j, x, *dftol, y);
00686     }
00687     //if (t_j == ldbl_INFTY)
00688     if (ldbl_EGlpNumIsEqqual (t_j, ldbl_INFTY))
00689       continue;
00690 
00691     if (ldbl_EGlpNumIsLess (t_j, t_max))
00692       ldbl_EGlpNumCopy (t_max, t_j);
00693   }
00694 
00695   if (ldbl_EGlpNumIsLeq (ldbl_INFTY, t_max))
00696   {
00697     rs->ratio_stat = RATIO_UNBOUNDED;
00698     ILL_CLEANUP;
00699   }
00700   /* if (t_max < 0.0) printf ("dIIhell\n"); */
00701 
00702   indx = -1;
00703   ldbl_EGlpNumZero (t_z);
00704   ldbl_EGlpNumZero (z_max);
00705   ldbl_EGlpNumZero (az_max);
00706 
00707   for (k = 0; k < lp->zA.nzcnt; k++)
00708   {
00709     zAj = &(lp->zA.coef[k]);
00710     ldbl_EGlpNumCopyAbs (azAj, *zAj);
00711     if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
00712       continue;
00713 
00714     ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
00715     j = lp->zA.indx[k];
00716     col = lp->nbaz[j];
00717 
00718     if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED)
00719       continue;
00720 
00721     ldbl_GET_XY_DRATIOTEST;
00722 
00723     if (ldbl_EGlpNumIsGreatZero (y) || lp->vstat[col] == STAT_ZERO)
00724       ldbl_EGlpNumCopyFrac (t_j, x, y);
00725 
00726     if (ldbl_EGlpNumIsLeq (t_j, t_max) && (ldbl_EGlpNumIsLess (az_max, azAj)))
00727     {
00728       ldbl_EGlpNumCopy (z_max, *zAj);
00729       ldbl_EGlpNumCopy (az_max, azAj);
00730       indx = j;
00731       ldbl_EGlpNumCopy (t_z, t_j);
00732     }
00733   }
00734 
00735 
00736   if (indx < 0)
00737   {
00738     rs->ratio_stat = RATIO_FAILED;
00739   }
00740   else
00741   {
00742     rs->eindex = indx;
00743     ldbl_EGlpNumCopy (rs->tz, t_z);
00744     ldbl_EGlpNumCopy (rs->pivotval, z_max);
00745     rs->ratio_stat = RATIO_BCHANGE;
00746 
00747     if (ldbl_EGlpNumIsLessZero (rs->tz))
00748     {
00749       ldbl_EGlpNumCopyAbs (rs->tz, t_max);
00750       ldbl_EGlpNumDivUiTo (rs->tz, 20);
00751       rs->coeffch = 1;
00752       ecol = lp->nbaz[indx];
00753       ldbl_EGlpNumCopyDiff (rs->ecoeff, lp->cz[ecol], lp->dz[indx]);
00754       switch (lp->vstat[ecol])
00755       {
00756       case STAT_LOWER:
00757         ldbl_EGlpNumAddInnProdTo (rs->ecoeff, rs->tz, az_max);
00758         break;
00759       case STAT_UPPER:
00760         ldbl_EGlpNumSubInnProdTo (rs->ecoeff, rs->tz, az_max);
00761         break;
00762       default:
00763         ldbl_EGlpNumZero (rs->tz);
00764         break;
00765       }
00766     }
00767   }
00768 
00769 CLEANUP:
00770   ldbl_ILLfct_update_counts (lp, CNT_DIIPIV, 0, rs->pivotval);
00771   ldbl_EGlpNumCopy (lp->upd.piv, rs->pivotval);
00772   ldbl_EGlpNumClearVar (x);
00773   ldbl_EGlpNumClearVar (y);
00774   ldbl_EGlpNumClearVar (t_j);
00775   ldbl_EGlpNumClearVar (z_max);
00776   ldbl_EGlpNumClearVar (t_max);
00777   ldbl_EGlpNumClearVar (t_z);
00778   ldbl_EGlpNumClearVar (az_max);
00779   ldbl_EGlpNumClearVar (azAj);
00780 }
00781 
00782 void ldbl_ILLratio_longdII_test (
00783   ldbl_lpinfo * lp,
00784   int lindex,
00785   int lvstat,
00786   ldbl_ratio_res * rs)
00787 {
00788   int j, k, indx = 0, tctr = 0;
00789   int col, ecol;
00790   int vs, bnd_exist = 0;
00791   int *perm = lp->upd.perm;
00792   int *ix = lp->upd.ix;
00793   int b_indx = -1;
00794   long double *t = lp->upd.t;
00795   long double *l,
00796     *u,
00797     *xb,
00798     *zAj = 0,
00799     x,
00800     y,
00801     t_j,
00802     z_max,
00803     t_max, t_z, theta, rcost, delta, zb_val, tb_val, az_max, azb_val, azAj;
00804   long double *pftol = &(lp->tol->pfeas_tol);
00805   long double *dftol = &(lp->tol->dfeas_tol);
00806   long double *pivtol = &(lp->tol->pivot_tol);
00807 
00808   ldbl_EGlpNumInitVar (x);
00809   ldbl_EGlpNumInitVar (azAj);
00810   ldbl_EGlpNumInitVar (y);
00811   ldbl_EGlpNumInitVar (t_j);
00812   ldbl_EGlpNumInitVar (z_max);
00813   ldbl_EGlpNumInitVar (az_max);
00814   ldbl_EGlpNumInitVar (t_max);
00815   ldbl_EGlpNumInitVar (t_z);
00816   ldbl_EGlpNumInitVar (theta);
00817   ldbl_EGlpNumInitVar (rcost);
00818   ldbl_EGlpNumInitVar (delta);
00819   ldbl_EGlpNumInitVar (zb_val);
00820   ldbl_EGlpNumInitVar (azb_val);
00821   ldbl_EGlpNumInitVar (tb_val);
00822   ldbl_EGlpNumZero (t_j);
00823   ldbl_EGlpNumZero (delta);
00824   ldbl_EGlpNumZero (zb_val);
00825   ldbl_EGlpNumZero (azb_val);
00826   ldbl_EGlpNumCopy (tb_val, ldbl_NINFTY);
00827 //#warning not sure about THIS line
00828   ldbl_EGlpNumZero (rs->pivotval);
00829 
00830   rs->coeffch = 0;
00831   rs->eindex = -1;
00832   rs->ratio_stat = RATIO_FAILED;
00833 
00834   ILL_IFTRACE2 ("%s:tctr %d\n", __func__, 0);
00835   lp->upd.tctr = 0;
00836   lp->upd.i = 0;
00837   ldbl_EGlpNumZero (lp->upd.tz);
00838   ldbl_EGlpNumZero (lp->upd.piv);
00839   ldbl_EGlpNumZero (lp->upd.c_obj);
00840   ldbl_EGlpNumZero (lp->upd.dty);
00841 
00842   xb = &(lp->xbz[lindex]);
00843   col = lp->baz[lindex];
00844   l = &(lp->lz[col]);
00845   u = &(lp->uz[col]);
00846   //rcost = (lvstat == STAT_LOWER) ? l - xb : xb - u;
00847   if (lvstat == STAT_LOWER)
00848     ldbl_EGlpNumCopyDiff (rcost, *l, *xb);
00849   else
00850     ldbl_EGlpNumCopyDiff (rcost, *xb, *u);
00851 
00852   for (k = 0, ldbl_EGlpNumCopy (t_max, ldbl_INFTY); k < lp->zA.nzcnt; k++)
00853   {
00854     zAj = &(lp->zA.coef[k]);
00855     if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
00856       continue;
00857 
00858     ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
00859     j = lp->zA.indx[k];
00860     col = lp->nbaz[j];
00861 
00862     if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED)
00863       continue;
00864     if (lp->vtype[col] == VBOUNDED)
00865     {
00866       bnd_exist++;
00867       continue;
00868     }
00869 
00870     ldbl_GET_XY_DRATIOTEST;
00871 
00872     if (ldbl_EGlpNumIsGreatZero (y))
00873     {
00874       //t_j = (x + dftol) / y;
00875 //#warning Using tolerances to add to result, is it right?
00876       ldbl_EGlpNumCopySum (t_j, x, *dftol);
00877       ldbl_EGlpNumDivTo (t_j, y);
00878     }
00879     else
00880     {
00881       if (lp->vstat[col] == STAT_ZERO)
00882         ldbl_EGlpNumCopyDiffRatio (t_j, x, *dftol, y);
00883     }
00884     if (ldbl_EGlpNumIsEqqual (t_j, ldbl_INFTY))
00885       continue;
00886 
00887     if (ldbl_EGlpNumIsLess (t_j, t_max))
00888       ldbl_EGlpNumCopy (t_max, t_j);
00889   }
00890   if (ldbl_EGlpNumIsLessZero (t_max))
00891   {
00892     /*printf ("dIIhell, %.4f\n", t_max); */
00893     rs->ratio_stat = RATIO_NEGATIVE;
00894     ILL_CLEANUP;
00895   }
00896 
00897   if (bnd_exist == 0 && ldbl_EGlpNumIsLeq (ldbl_INFTY, t_max))
00898   {
00899     rs->ratio_stat = RATIO_UNBOUNDED;
00900     /*
00901      * printf ("x = %.8f, b = %.2f \n", lp->xbz[lindex], (lvstat == STAT_LOWER ) ? lp->lz[lp->baz[lindex]] : lp->uz[lp->baz[lindex]]);
00902      */
00903     ILL_CLEANUP;
00904   }
00905 
00906   if (bnd_exist != 0)
00907   {
00908     for (k = 0; k < lp->zA.nzcnt; k++)
00909     {
00910       zAj = &(lp->zA.coef[k]);
00911       if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
00912         continue;
00913 
00914       ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
00915       j = lp->zA.indx[k];
00916       col = lp->nbaz[j];
00917 
00918       if (lp->vtype[col] != VBOUNDED)
00919         continue;
00920 
00921       ldbl_GET_XY_DRATIOTEST;
00922 
00923       if (ldbl_EGlpNumIsGreatZero (y))
00924       {
00925         ldbl_EGlpNumCopyFrac (t_j, x, y);
00926         if (ldbl_EGlpNumIsLeq (t_j, t_max))
00927         {
00928           ldbl_EGlpNumCopy (t[tctr], t_j);
00929           ix[tctr] = k;
00930           tctr++;
00931         }
00932       }
00933     }
00934   }
00935 
00936   if (tctr != 0)
00937   {
00938     for (j = 0; j < tctr; j++)
00939       perm[j] = j;
00940     ldbl_ILLutil_EGlpNum_perm_quicksort (perm, t, tctr);
00941 
00942     for (j = 0; j < tctr; j++)
00943     {
00944 
00945       ldbl_EGlpNumCopy (t_j, t[perm[j]]);
00946       /* we use x as temporal storage */
00947       //lp->upd.c_obj += (t_j - delta) * rcost;
00948       ldbl_EGlpNumCopy (x, t_j);
00949       ldbl_EGlpNumSubTo (x, delta);
00950       ldbl_EGlpNumAddInnProdTo (lp->upd.c_obj, x, rcost);
00951       ldbl_EGlpNumCopy (delta, t_j);
00952        /*HHH*/ k = ix[perm[j]];
00953       zAj = &(lp->zA.coef[k]);
00954       indx = lp->zA.indx[k];
00955       col = lp->nbaz[indx];
00956       l = &(lp->lz[col]);
00957       u = &(lp->uz[col]);
00958       vs = lp->vstat[col];
00959       //theta = (vs == STAT_UPPER) ? (l - u) * zAj : (u - l) * zAj;
00960       ldbl_EGlpNumCopyDiff (theta, *l, *u);
00961       ldbl_EGlpNumMultTo (theta, *zAj);
00962       if (vs != STAT_UPPER)
00963         ldbl_EGlpNumSign (theta);
00964       if (lvstat == STAT_LOWER)
00965         ldbl_EGlpNumAddTo (rcost, theta);
00966       else
00967         ldbl_EGlpNumSubTo (rcost, theta);
00968 
00969       if (ldbl_EGlpNumIsLeq (rcost, *pftol))
00970       {
00971         rs->eindex = indx;
00972         ldbl_EGlpNumCopy (rs->tz, t_j);
00973         ldbl_EGlpNumCopy (rs->pivotval, *zAj);
00974         rs->ratio_stat = RATIO_BCHANGE;
00975 
00976         if (ldbl_EGlpNumIsLessZero (rs->tz))
00977         {
00978           ldbl_EGlpNumZero (rs->tz);
00979           rs->coeffch = 1;
00980           //rs->ecoeff = lp->cz[col] - lp->dz[indx];
00981           ldbl_EGlpNumCopyDiff (rs->ecoeff, lp->cz[col], lp->dz[indx]);
00982           //lp->upd.c_obj += (rs->tz - delta) * rcost; note ts->tz == 0;
00983           ldbl_EGlpNumSubInnProdTo (lp->upd.c_obj, delta, rcost);
00984         }
00985         ILL_IFTRACE2 ("%s:tctr %d\n", __func__, tctr);
00986         lp->upd.tctr = tctr;
00987         lp->upd.i = j;
00988         ldbl_EGlpNumCopy (lp->upd.tz, rs->tz);
00989         ILL_CLEANUP;
00990       }
00991     }
00992     ILL_IFTRACE2 ("%s:tctr %d\n", __func__, tctr);
00993     lp->upd.tctr = tctr;
00994     lp->upd.i = tctr;
00995     ldbl_EGlpNumCopy (lp->upd.tz, t_j);
00996     ldbl_EGlpNumCopy (zb_val, *zAj);
00997     ldbl_EGlpNumCopyAbs (azb_val, zb_val);
00998     ldbl_EGlpNumCopy (tb_val, t_j);
00999     b_indx = indx;
01000   }
01001 
01002   if (bnd_exist != 0 && ldbl_EGlpNumIsLeq (ldbl_INFTY, t_max))
01003   {
01004     rs->ratio_stat = RATIO_UNBOUNDED;
01005     /* printf ("rcost: %.8f\n", rcost); */
01006     ILL_CLEANUP;
01007   }
01008 
01009   ldbl_EGlpNumZero (z_max);
01010   ldbl_EGlpNumZero (az_max);
01011   indx = -1;
01012   ldbl_EGlpNumZero (t_z);
01013   for (k = 0; k < lp->zA.nzcnt; k++)
01014   {
01015     zAj = &(lp->zA.coef[k]);
01016     ldbl_EGlpNumCopyAbs (azAj, *zAj);
01017     if (!ldbl_EGlpNumIsNeqZero (*zAj, *pivtol))
01018       continue;
01019 
01020     ldbl_EGlpNumCopy (t_j, ldbl_INFTY);
01021     j = lp->zA.indx[k];
01022     col = lp->nbaz[j];
01023 
01024     if (lp->vtype[col] == VARTIFICIAL || lp->vtype[col] == VFIXED ||
01025         lp->vtype[col] == VBOUNDED)
01026       continue;
01027 
01028     ldbl_GET_XY_DRATIOTEST;
01029 
01030     if (ldbl_EGlpNumIsGreatZero (y) || lp->vstat[col] == STAT_ZERO)
01031       ldbl_EGlpNumCopyFrac (t_j, x, y);
01032 
01033     if (ldbl_EGlpNumIsLeq (t_j, t_max))
01034     {
01035       if (ldbl_EGlpNumIsLess (az_max, azAj))
01036       {
01037         ldbl_EGlpNumCopy (z_max, *zAj);
01038         ldbl_EGlpNumCopy (az_max, azAj);
01039         indx = j;
01040         ldbl_EGlpNumCopy (t_z, t_j);
01041       }
01042     }
01043   }
01044 
01045   if (indx < 0)
01046   {
01047     rs->ratio_stat = RATIO_FAILED;
01048     ILL_CLEANUP;
01049   }
01050   if ((tctr == 0) || (ldbl_EGlpNumIsLessZero (tb_val)) ||
01051       (tctr != 0 && ldbl_EGlpNumIsLeq (tb_val, t_z) &&
01052        ldbl_EGlpNumIsLeq (azb_val, az_max)))
01053   {
01054     /* we use x as temporal vvariable */
01055     /* lp->upd.c_obj += (t_z - delta) * rcost; */
01056     ldbl_EGlpNumCopyDiff (x, t_z, delta);
01057     ldbl_EGlpNumAddInnProdTo (lp->upd.c_obj, x, rcost);
01058     ldbl_EGlpNumCopy (delta, t_z);
01059     rs->eindex = indx;
01060     ldbl_EGlpNumCopy (rs->tz, t_z);
01061     ldbl_EGlpNumCopy (rs->pivotval, z_max);
01062     rs->ratio_stat = RATIO_BCHANGE;
01063   }
01064   /* For now */
01065   else if (tctr != 0)
01066   {
01067     rs->eindex = b_indx;
01068     ldbl_EGlpNumCopy (rs->tz, tb_val);
01069     ldbl_EGlpNumCopy (rs->pivotval, zb_val);
01070     rs->ratio_stat = RATIO_BCHANGE;
01071     lp->upd.i -= 1;
01072   }
01073 
01074   if (ldbl_EGlpNumIsLessZero (rs->tz))
01075   {
01076     /* if (tctr != 0) printf ("despite long step\n"); */
01077     /* rs->tz = fabs (t_max / 20.0); */
01078     ldbl_EGlpNumCopyAbs (rs->tz, t_max);
01079     ldbl_EGlpNumDivUiTo (rs->tz, 20);
01080     rs->coeffch = 1;
01081 
01082     ecol = lp->nbaz[indx];
01083     if (lp->vstat[ecol] == STAT_LOWER)
01084     {
01085       /*rs->ecoeff = lp->cz[ecol] - lp->dz[indx] + rs->tz * fabs (z_max); */
01086       ldbl_EGlpNumCopy (rs->ecoeff, az_max);
01087       ldbl_EGlpNumMultTo (rs->ecoeff, rs->tz);
01088       ldbl_EGlpNumAddTo (rs->ecoeff, lp->cz[ecol]);
01089       ldbl_EGlpNumSubTo (rs->ecoeff, lp->dz[indx]);
01090     }
01091     else if (lp->vstat[ecol] == STAT_UPPER)
01092     {
01093       /*rs->ecoeff = lp->cz[ecol] - lp->dz[indx] - rs->tz * fabs (z_max); */
01094       ldbl_EGlpNumCopy (rs->ecoeff, az_max);
01095       ldbl_EGlpNumMultTo (rs->ecoeff, rs->tz);
01096       ldbl_EGlpNumSign (rs->ecoeff);
01097       ldbl_EGlpNumAddTo (rs->ecoeff, lp->cz[ecol]);
01098       ldbl_EGlpNumSubTo (rs->ecoeff, lp->dz[indx]);
01099     }
01100     else
01101     {
01102       /*rs->ecoeff = lp->cz[ecol] - lp->dz[indx]; */
01103       ldbl_EGlpNumCopyDiff (rs->ecoeff, lp->cz[ecol], lp->dz[indx]);
01104       ldbl_EGlpNumZero (rs->tz);
01105     }
01106     /* we use x as temporal storage */
01107     /*lp->upd.c_obj += (rs->tz - delta) * rcost; */
01108     ldbl_EGlpNumCopy (x, rs->tz);
01109     ldbl_EGlpNumSubTo (x, delta);
01110     ldbl_EGlpNumAddInnProdTo (lp->upd.c_obj, x, rcost);
01111   }
01112 
01113 CLEANUP:
01114   ldbl_ILLfct_update_counts (lp, CNT_DIIPIV, 0, rs->pivotval);
01115   ldbl_EGlpNumCopy (lp->upd.piv, rs->pivotval);
01116   ldbl_EGlpNumClearVar (x);
01117   ldbl_EGlpNumClearVar (y);
01118   ldbl_EGlpNumClearVar (t_j);
01119   ldbl_EGlpNumClearVar (z_max);
01120   ldbl_EGlpNumClearVar (az_max);
01121   ldbl_EGlpNumClearVar (t_max);
01122   ldbl_EGlpNumClearVar (t_z);
01123   ldbl_EGlpNumClearVar (theta);
01124   ldbl_EGlpNumClearVar (rcost);
01125   ldbl_EGlpNumClearVar (delta);
01126   ldbl_EGlpNumClearVar (zb_val);
01127   ldbl_EGlpNumClearVar (azb_val);
01128   ldbl_EGlpNumClearVar (tb_val);
01129   ldbl_EGlpNumClearVar (azAj);
01130 }
01131 
01132 void ldbl_ILLratio_pivotin_test (
01133   ldbl_lpinfo * lp,
01134   int *rlist,
01135   int rcnt,
01136   ldbl_ratio_res * rs)
01137 {
01138   int i, k, col;
01139   long double *x, *l, *u;
01140   long double ay_ij,
01141     at_i, at_l, at_u, ayi_max, y_ij, t_i, t_l, t_u, t_max, yi_max;
01142   long double *pivtol = &(lp->tol->pivot_tol);
01143 
01144   if (rcnt <= 0 || rs == NULL)
01145     return;
01146   ldbl_EGlpNumInitVar (ay_ij);
01147   ldbl_EGlpNumInitVar (at_i);
01148   ldbl_EGlpNumInitVar (at_l);
01149   ldbl_EGlpNumInitVar (at_u);
01150   ldbl_EGlpNumInitVar (ayi_max);
01151   ldbl_EGlpNumInitVar (t_max);
01152   ldbl_EGlpNumInitVar (y_ij);
01153   ldbl_EGlpNumInitVar (t_i);
01154   ldbl_EGlpNumInitVar (t_l);
01155   ldbl_EGlpNumInitVar (t_u);
01156   ldbl_EGlpNumInitVar (yi_max);
01157   rs->boundch = 0;
01158   rs->lindex = -1;
01159   ldbl_EGlpNumZero (rs->tz);
01160   rs->ratio_stat = RATIO_FAILED;
01161   rs->lvstat = -1;
01162   ldbl_EGlpNumZero (rs->pivotval);
01163   ldbl_EGlpNumZero (rs->lbound);
01164 
01165   for (i = 0; i < rcnt; i++)
01166     lp->iwork[rlist[i]] = 1;
01167 
01168   for (k = 0, ldbl_EGlpNumCopy (t_max, ldbl_INFTY); k < lp->yjz.nzcnt; k++)
01169   {
01170     ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
01171     if (!ldbl_EGlpNumIsNeqZero (y_ij, *pivtol))
01172       continue;
01173 
01174     i = lp->yjz.indx[k];
01175     if (lp->iwork[lp->baz[i]] == 1)
01176       continue;
01177     x = &(lp->xbz[i]);
01178     col = lp->baz[i];
01179     l = &(lp->lz[col]);
01180     u = &(lp->uz[col]);
01181     ldbl_EGlpNumCopy (t_u, ldbl_INFTY);
01182     ldbl_EGlpNumCopy (at_u, ldbl_INFTY);
01183     ldbl_EGlpNumCopy (t_l, ldbl_NINFTY);
01184     ldbl_EGlpNumCopy (at_l, ldbl_INFTY);
01185 
01186     if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
01187     {
01188       ldbl_EGlpNumCopyDiffRatio (t_l, *x, *l, y_ij);
01189       ldbl_EGlpNumCopyAbs (at_l, t_l);
01190       if (ldbl_EGlpNumIsLess (at_l, t_max))
01191         ldbl_EGlpNumCopy (t_max, at_l);
01192     }
01193     if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
01194     {
01195       ldbl_EGlpNumCopyDiffRatio (t_u, *x, *u, y_ij);
01196       ldbl_EGlpNumCopyAbs (at_u, t_u);
01197       if (ldbl_EGlpNumIsLess (at_u, t_max))
01198         ldbl_EGlpNumCopy (t_max, at_u);
01199     }
01200   }
01201 
01202   if (ldbl_EGlpNumIsLeq (ldbl_INFTY, t_max))
01203   {
01204     rs->ratio_stat = RATIO_UNBOUNDED;
01205     ILL_CLEANUP;
01206   }
01207 
01208   ldbl_EGlpNumZero (yi_max);
01209   ldbl_EGlpNumZero (ayi_max);
01210   ldbl_EGlpNumMultUiTo (t_max, 101);
01211   ldbl_EGlpNumDivUiTo (t_max, 100);
01212   for (k = 0; k < lp->yjz.nzcnt; k++)
01213   {
01214     ldbl_EGlpNumCopy (y_ij, lp->yjz.coef[k]);
01215     ldbl_EGlpNumCopyAbs (ay_ij, y_ij);
01216     if (!ldbl_EGlpNumIsNeqZero (y_ij, *pivtol))
01217       continue;
01218 
01219     i = lp->yjz.indx[k];
01220     if (lp->iwork[lp->baz[i]] == 1)
01221       continue;
01222     x = &(lp->xbz[i]);
01223     col = lp->baz[i];
01224     l = &(lp->lz[col]);
01225     u = &(lp->uz[col]);
01226 
01227     ldbl_EGlpNumCopy (t_u, ldbl_INFTY);
01228     ldbl_EGlpNumCopy (at_u, t_u);
01229     ldbl_EGlpNumCopy (t_l, ldbl_NINFTY);
01230     ldbl_EGlpNumCopy (at_l, t_u);
01231     if (ldbl_EGlpNumIsNeqq (*l, ldbl_NINFTY))
01232     {
01233       ldbl_EGlpNumCopyDiffRatio (t_l, *x, *l, y_ij);
01234       ldbl_EGlpNumCopyAbs (at_l, t_l);
01235     }
01236     if (ldbl_EGlpNumIsNeqq (*u, ldbl_INFTY))
01237     {
01238       ldbl_EGlpNumCopyDiffRatio (t_u, *x, *u, y_ij);
01239       ldbl_EGlpNumCopyAbs (at_u, t_u);
01240     }
01241     //t_i = (fabs (t_l) < fabs (t_u)) ? t_l : t_u;
01242     if (ldbl_EGlpNumIsLess (at_l, at_u))
01243     {
01244       ldbl_EGlpNumCopy (t_i, t_l);
01245       ldbl_EGlpNumCopy (at_i, at_l);
01246     }
01247     else
01248     {
01249       ldbl_EGlpNumCopy (t_i, t_u);
01250       ldbl_EGlpNumCopy (at_i, at_u);
01251     }
01252     /*if (fabs (t_i) <= t_max + t_max * (1.0e-2)) */
01253     if (ldbl_EGlpNumIsLeq (at_i, t_max))
01254     {
01255       if (ldbl_EGlpNumIsLess (ayi_max, ay_ij))
01256       {
01257         ldbl_EGlpNumCopy (yi_max, y_ij);
01258         ldbl_EGlpNumCopy (ayi_max, ay_ij);
01259         rs->lindex = i;
01260         ldbl_EGlpNumCopy (rs->tz, t_i);
01261         rs->lvstat = (ldbl_EGlpNumIsLess (at_l, at_u)) ? STAT_LOWER : STAT_UPPER;
01262       }
01263     }
01264   }
01265 
01266   if (rs->lindex < 0)
01267   {
01268     rs->ratio_stat = RATIO_FAILED;
01269   }
01270   else
01271   {
01272     rs->ratio_stat = RATIO_BCHANGE;
01273     ldbl_EGlpNumCopy (rs->pivotval, yi_max);
01274   }
01275 CLEANUP:
01276   for (i = 0; i < rcnt; i++)
01277     lp->iwork[rlist[i]] = 0;
01278   ldbl_EGlpNumClearVar (t_max);
01279   ldbl_EGlpNumClearVar (ay_ij);
01280   ldbl_EGlpNumClearVar (at_i);
01281   ldbl_EGlpNumClearVar (at_l);
01282   ldbl_EGlpNumClearVar (at_u);
01283   ldbl_EGlpNumClearVar (ayi_max);
01284   ldbl_EGlpNumClearVar (y_ij);
01285   ldbl_EGlpNumClearVar (t_i);
01286   ldbl_EGlpNumClearVar (t_l);
01287   ldbl_EGlpNumClearVar (t_u);
01288   ldbl_EGlpNumClearVar (yi_max);
01289   return;
01290 }

Generated on Thu Mar 29 09:32:31 2012 for QSopt_ex by  doxygen 1.4.7