2 #ifndef COIN_OSL_C_INCLUDE
9 #define COIN_OSL_C_INCLUDE
14 #define C_EKK_GO_SPARSE 200
18 #if __BYTE_ORDER == __LITTLE_ENDIAN
39 int *mpt,
int first_nonzero);
42 int *mpt,
int ipivrw,
int *spare);
46 int *mpt2,
double dalpha,
int orig_nincol,
47 int npivot,
int *nuspikp,
48 const int ipivrw,
int *spare);
51 double *dpermu,
int *mpt,
int numberNonZero);
54 double *dwork1,
int *mpt,
int *nincolp);
56 double *dpermu1,
int *mpt1,
int *nincolp,
57 double *dwork1_ft,
int *mpt_ft,
int *nincolp_ft);
61 inline void c_ekkscpy(
int n,
const int *marr1,
int *marr2)
65 inline void c_ekkdcpy(
int n,
const double *marr1,
double *marr2)
73 void c_ekkzero(
int length,
int n,
void *array);
90 #define c_ekkscpy_0_1(s, ival, array) CoinFillN(array, s, ival)
91 #define c_ekks1cpy(n, marr1, marr2) CoinMemcpyN(marr1, n, marr2)
95 int *
clp_int(
int number_entries);
99 #define SLACK_VALUE -1.0
100 #define C_EKK_REMOVE_LINK(hpiv, hin, link, ipivot) \
102 int ipre = link[ipivot].pre; \
103 int isuc = link[ipivot].suc; \
105 link[ipre].suc = isuc; \
108 hpiv[hin[ipivot]] = isuc; \
111 link[isuc].pre = ipre; \
115 #define C_EKK_ADD_LINK(hpiv, nzi, link, npr) \
117 int ifiri = hpiv[nzi]; \
119 link[npr].suc = ifiri; \
122 link[ifiri].pre = npr; \
128 #define SHIFT_INDEX(limit) (limit)
129 #define UNSHIFT_INDEX(limit) (limit)
130 #define SHIFT_REF(arr, ind) (arr)[ind]
134 #define SHIFT_INDEX(limit) ((limit) << 3)
135 #define UNSHIFT_INDEX(limit) ((unsigned int)(limit) >> 3)
136 #define SHIFT_REF(arr, ind) (*(double *)((char *)(arr) + (ind)))
141 #define NOT_ZERO(x) (((*((reinterpret_cast< unsigned char * >(&x)) + 7)) & 0x7F) != 0)
143 #define NOT_ZERO(x) ((x) != 0.0)
146 #define SWAP(type, _x, _y) \
153 #define UNROLL_LOOP_BODY1(code) \
159 #define UNROLL_LOOP_BODY2(code) \
165 #define UNROLL_LOOP_BODY4(code) \
167 { code } { code } { code } { \
180 int *hcoli = fact->xecadr;
181 double *dluval = fact->xeeadr;
182 double *dvalpv = fact->kw3adr;
183 int *mrstrt = fact->xrsadr;
184 int *hrowi = fact->xeradr;
185 int *mcstrt = fact->xcsadr;
186 int *hinrow = fact->xrnadr;
187 int *hincol = fact->xcnadr;
188 int *hpivro = fact->krpadr;
189 int *hpivco = fact->kcpadr;
191 int nnentl = fact->nnentl;
192 int nnentu = fact->nnentu;
193 int kmxeta = fact->kmxeta;
194 int xnewro = *xnewrop;
195 int ncompactions = *ncompactionsp;
197 MACTION_T *maction =
reinterpret_cast< MACTION_T *
>(maction_void);
202 int jj, kk, kr, nz, jj1, jj2, kce, kcs, kqq, npr;
206 int knpr, irow, iadd32, ibase;
212 int kipie, kcpiv, knprs, knpre;
214 double multip, elemnt;
215 int ipivot, jpivot, epivro, epivco, lstart, nfirst;
216 int nzpivj, kfill, kstart;
220 int noRoomForDense = 0;
221 int if_sparse_update = fact->if_sparse_update;
225 const int nrow = fact->nrow;
231 lstart = nnetas - nnentl + 1;
232 for (i = lstart; i <= nnetas; ++i) {
236 for (i = 1; i <= nrow; ++i) {
238 mwork[i].pre = i - 1;
239 mwork[i].suc = i + 1;
248 for (count = 1; count <= nrow; ++count) {
251 if (!(hpivco[1] <= 0)) {
252 int small_pivot = c_ekkcsin(fact,
258 if (fact->invok >= 0) {
262 if (fact->npivots >= nrow) {
268 if (!(hpivro[1] <= 0)) {
269 irtcod = c_ekkrsin(fact,
276 &kmxeta, &ncompactions,
279 if (irtcod < 0 || fact->invok >= 0) {
286 if (fact->npivots >= nrow) {
289 lstart = nnetas - nnentl + 1;
293 irtcod = c_ekkfpvt(fact,
295 nsingp, xrejctp, &ipivot, &jpivot);
301 c_ekkprpv(fact, rlink, clink,
302 *xrejctp, ipivot, jpivot);
304 epivco = hincol[jpivot];
306 mcstrt[fact->xnetal] = lstart - 1;
307 hpivco[fact->xnetal] = ipivot;
308 epivro = hinrow[ipivot];
310 kipis = mrstrt[ipivot];
311 pivot = dluval[kipis];
313 kipie = kipis + epivr1;
317 double size = nrow - fact->npivots;
318 if (size > GO_DENSE && (nnentu - fact->nuspike) * GO_DENSE_RATIO > size * size) {
329 for (k = kipis; k <= kipie; ++k) {
331 dvalpv[k - kipis + 1] = dluval[k];
332 maction[irow] =
static_cast< MACTION_T
>(k - kipis + 1);
336 kcpiv = mcstrt[jpivot] - 1;
337 for (nzpivj = 1; nzpivj <= epivco; ++nzpivj) {
350 knpre = knprs + enpr;
356 if (jpivot == hcoli[knpr]) {
363 multip = -dluval[knpr] * dpivx;
366 dluval[knpr] = dluval[knpre];
367 hcoli[knpr] = hcoli[knpre];
372 kfill = epivr1 - (knpre - knprs + 1);
373 nres = ((knpre - knprs + 1) & 1) + knprs;
380 if (maction[j] == 0) {
384 maction[j] =
static_cast< MACTION_T
>(-maction[j]);
385 dluval[knprs] += multip * dvalpv[jj];
386 d1 = fabs(dluval[knprs]);
389 j2 = hcoli[nres + 1];
391 for (kr = nres; kr < knpre; kr += 2) {
396 maction[j1] =
static_cast< MACTION_T
>(-maction[j1]);
397 dluval[kr] += multip * dvalpv[jj1];
398 cancel = cancel || !(fact->zeroTolerance < d1);
399 d1 = fabs(dluval[kr]);
405 maction[j2] =
static_cast< MACTION_T
>(-maction[j2]);
406 dluval[kr + 1] += multip * dvalpv[jj2];
407 cancel = cancel || !(fact->zeroTolerance < d1);
408 d1 = fabs(dluval[kr + 1]);
413 cancel = cancel || !(fact->zeroTolerance < d1);
430 kfill = epivr1 - (knpre - knprs + 1);
433 for (kr = knprs; kr <= knpre; kr++) {
441 maction[j1] = -maction[j1];
442 dluval[kr] += multip * dvalpv[jj1];
443 d1 = fabs(dluval[kr]);
444 cancel = cancel || !(fact->zeroTolerance < d1);
454 for (kr = knprs; kr <= knpre; ++kr) {
456 if (fabs(dluval[kr]) > fact->zeroTolerance) {
458 dluval[kstart] = dluval[kr];
466 kce = kcs + hincol[j];
467 for (kk = kcs; kk <= kce; ++kk) {
468 if (hrowi[kk] == npr) {
469 hrowi[kk] = hrowi[kce];
481 for (k = kipis; k <= kipie; ++k) {
482 maction[hcoli[k]] =
static_cast< MACTION_T
>(-maction[hcoli[k]]);
485 naft = mwork[npr].suc;
486 kqq = mrstrt[naft] - knpre - 1;
492 if (!(xnewro + nznpr + 1 < lstart)) {
493 if (!(nnentu + nznpr + 1 < lstart)) {
503 int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
504 kmxeta += xnewro - iput;
509 kipis = mrstrt[ipivot] + 1;
510 kipie = kipis + epivr1 - 1;
517 knpre = knprs + enpr - 1;
525 npre = mwork[npr].pre;
530 mwork[naft].pre = npre;
531 mwork[npre].suc = naft;
533 mwork[nfirst].pre = npr;
534 mwork[nlast].suc = npr;
535 mwork[npr].pre = nlast;
536 mwork[npr].suc = nfirst;
539 mrstrt[npr] = kstart + 1;
540 nmove = knpre - knprs + 1;
541 ibase = kstart + 1 - knprs;
542 for (kr = knprs; kr <= knpre; ++kr) {
543 dluval[ibase + kr] = dluval[kr];
544 hcoli[ibase + kr] = hcoli[kr];
565 if (kstart + (nrow << 1) + 100 < lstart) {
566 ileft = ((nrow - fact->npivots + 32) & -32);
567 if (kstart + ileft * ileft + 32 < lstart) {
569 xnewro =
CoinMax(kstart, xnewro);
570 xnewro = (xnewro & -32) + ileft;
572 xnewro = ((kstart + 31) & -32);
580 }
else if (!(nnentu + kqq + 2 < lstart)) {
585 for (kr = kipis; kr <= kipie; ++kr) {
589 elemnt = multip * dvalpv[jj];
590 if (fabs(elemnt) > fact->zeroTolerance) {
592 dluval[kstart] = elemnt;
601 if (xnewco + 1 >= lstart) {
602 if (xnewco + nz + 1 >= lstart) {
604 if (nnentu + nz + 1 < lstart) {
605 xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
608 kcpiv = mcstrt[jpivot] - 1;
619 mcstrt[j] = xnewco + 1;
620 ibase = mcstrt[j] - kcs;
621 for (kk = kcs; kk <= kce; ++kk) {
622 hrowi[ibase + kk] = hrowi[kk];
625 kce = xnewco + kce - kcs + 1;
630 }
else if (hrowi[kce + 1] != 0) {
632 if (xnewco + nz + 1 >= lstart) {
634 if (nnentu + nz + 1 < lstart) {
635 xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
638 kcpiv = mcstrt[jpivot] - 1;
649 mcstrt[j] = xnewco + 1;
650 ibase = mcstrt[j] - kcs;
651 for (kk = kcs; kk <= kce; ++kk) {
652 hrowi[ibase + kk] = hrowi[kk];
655 kce = xnewco + kce - kcs + 1;
659 hrowi[kce + 1] = npr;
663 maction[j] =
static_cast< MACTION_T
>(-maction[j]);
670 hinrow[npr] = kstart - mrstrt[npr] + 1;
672 if (!(xnewco + 1 < lstart)) {
673 xnewco = c_ekkclco(fact, hrowi, mcstrt, hincol, xnewco);
676 kcpiv = mcstrt[jpivot] - 1;
678 if (!(xnewro + 1 < lstart)) {
679 int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
680 kmxeta += xnewro - iput;
684 kipis = mrstrt[ipivot] + 1;
685 kipie = kipis + epivr1 - 1;
691 dluval[lstart] = multip;
694 #define INLINE_AFPV 3
715 double d = fabs(els[j]);
725 double maxaij = fabs(els[0]);
726 for (j = 1; j < nel; j++) {
727 double d = fabs(els[j]);
738 if ((nel & 1) != 0) {
739 maxaij = fabs(els[0]);
745 double d = fabs(els[j]);
754 SWAP(
int, index[koff], index[0]);
755 SWAP(
double, els[koff], els[0]);
760 int nzi = hinrow[npr];
769 int nn = mcstrt[fact->xnetal] - lstart + 1;
770 c_ekkafpv(hrowi + lstart, hcoli, dluval, mrstrt, hinrow, nn);
774 for (k = kipis; k <= kipie; ++k) {
775 maction[hcoli[k]] = 0;
779 for (k = kipis; k <= kipie; ++k) {
782 if (!(nzj <= 0) && !((clink[j].pre > nrow && nzj != 1))) {
787 for (k = kipis; k <= kipie; ++k) {
795 fact->nuspike += hinrow[ipivot];
800 int ndense = nrow - fact->npivots;
801 if (!(xnewro + ndense * ndense >= lstart)) {
804 c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
806 for (i = 1; i <= nrow; ++i) {
807 if (clink[i].pre >= 0) {
809 maction[i] =
static_cast< short int >(iput);
814 for (i = 1; i <= nrow; ++i) {
815 if (rlink[i].pre >= 0) {
816 nspare = nspare + ndense - hinrow[i];
819 if (iput != nrow - fact->npivots) {
821 c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
824 int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
825 kmxeta += xnewro - iput;
830 if (xnewro + nspare + ndense * ndense >= lstart) {
831 c_ekkizero(nrow, reinterpret_cast< int * >(maction + 1));
834 c_ekkrwct(fact, dluval, hcoli, mrstrt, hinrow, mwork,
835 rlink, maction, dvalpv,
838 if (nnentu + nnentl > nrow * 5 && (ndense * ndense) > (nnentu + nnentl) >> 2 && !if_sparse_update) {
839 fact->ndenuc = ndense;
841 irtcod = c_ekkcmfd(fact,
842 (reinterpret_cast< int * >(dvalpv) + 1),
844 (reinterpret_cast< int * >(maction + 1)) + 1,
858 if (!noRoomForDense) {
859 int etasize =
CoinMax(4 * fact->nnentu + (nnetas - fact->nnentl) + 1000, fact->eta_size);
860 noRoomForDense = ndense;
861 fact->eta_size =
CoinMin(static_cast< int >(1.2 * fact->eta_size), etasize);
862 if (fact->maxNNetas > 0 && fact->eta_size > fact->maxNNetas) {
863 fact->eta_size = fact->maxNNetas;
873 int iput = c_ekkrwcs(fact, dluval, hcoli, mrstrt, hinrow, mwork, nfirst);
874 kmxeta += xnewro - iput;
881 mwork[nrow + 1].pre = nfirst;
882 mwork[nrow + 1].suc = nlast;
884 fact->nnentl = nnentl;
885 fact->nnentu = nnentu;
886 fact->kmxeta = kmxeta;
888 *ncompactionsp = ncompactions;
void clp_free(void *oldArray)
void c_ekkscpy(int n, const int *marr1, int *marr2)
#define SWAP(type, _x, _y)
int * clp_int(int number_entries)
int c_ekketsj(registerEKKfactinfo *fact, double *dwork1, int *mpt2, double dalpha, int orig_nincol, int npivot, int *nuspikp, const int ipivrw, int *spare)
#define UNROLL_LOOP_BODY2(code)
void c_ekk_Set(int *array, int bit)
void clp_setup_pointers(EKKfactinfo *fact)
int c_ekkbtrn(register const EKKfactinfo *fact, double *dwork1, int *mpt, int first_nonzero)
int c_ekkftrn(register const EKKfactinfo *fact, double *dwork1, double *dpermu, int *mpt, int numberNonZero)
void CoinMemcpyN(const T *from, const CoinBigIndex size, T *to)
This helper function copies an array to another location.
double * clp_double(int number_entries)
int c_ekkbtrn_ipivrw(register const EKKfactinfo *fact, double *dwork1, int *mpt, int ipivrw, int *spare)
#define SHIFT_INDEX(limit)
void * clp_malloc(int number_entries)
void c_ekk_Unset(int *array, int bit)
int c_ekklfct(register EKKfactinfo *fact)
void c_ekkzero(int length, int n, void *array)
void c_ekkczero(int n, char *marray)
void CoinZeroN(T *to, const CoinBigIndex size)
This helper function fills an array with zero.
T CoinMax(const T x1, const T x2)
Return the larger (according to operator<() of the arguments.
void c_ekkdcpy(int n, const double *marr1, double *marr2)
void c_ekkftrn2(register EKKfactinfo *fact, double *dwork1, double *dpermu1, int *mpt1, int *nincolp, double *dwork1_ft, int *mpt_ft, int *nincolp_ft)
int c_ekkftrn_ft(register EKKfactinfo *fact, double *dwork1, int *mpt, int *nincolp)
void c_ekkizero(int n, int *marray)
void clp_memory(int type)
T CoinMin(const T x1, const T x2)
Return the smaller (according to operator<() of the arguments.
#define UNROLL_LOOP_BODY4(code)
int c_ekk_IsSet(const int *array, int bit)
int c_ekkslcf(register const EKKfactinfo *fact)
#define C_EKK_ADD_LINK(hpiv, nzi, link, npr)
void c_ekkdzero(int n, double *marray)