Re: Patches to improve pseudocost initialisatiion speed

Chris Matrakidis <[email protected]>
Newsgroups gmane.comp.gnu.glpk
Message-ID <CAEaJ3S73+nYrpeGTnfz8QYW-=yZZg8fdx+Od35cGpHXcoHyStA@mail.gmail.com>
Hi Andrew,

Leaving aside the issue of the API for the time being, I did some
additional experiments in this area.

> To resolve this
> issue some time ago I implemented another factorization of the basis
> based on Schur complement (please see comments in src/bflib/scf.h),
> where B0 = L0 * U0 is not changed, so B = B0 can be easily restored as
> in case of "eta file".

I did try this, and it seems to work OK. The attached two patches
implement it on top of the three original patches I sent (pcost1.patch
to pcost3.patch). The first one just makes the pseudocost
initialisation use Schur complement updates, to use as a baseline
since the pseudocosts calculated may be different now. The second ones
introduces a bfd_reset() function that restores the original
factorisation and uses this instead of bfd_copy().

Unfortunately, it looks like that the speed gain from using bfd_reset
does not offset the overhead of Schur complement updates, so
Forrest-Tomlin update with bfd_copy seems to be the preferred option
here.

> This approach, however, limits
> the number of simplex iteration (say, to 100-200), since if
> refactorization is needed, to keep B0 = L0 * U0 we should stop.

I suppose this means that a new simplex flag may be needed, to stop
when re-factorisation is required.

Best Regards,

Chris Matrakidis

_______________________________________________
Help-glpk mailing list
[email protected]
https://lists.gnu.org/mailman/listinfo/help-glpk
schur1.patch (text/x-patch, 1 KB)
commit 9734f4815fda152a20491a617818badb631adbf6
Author: Chris Matrakidis <[email protected]>
Date:   Sat Jan 14 16:24:45 2017 +0200

    use schur complement update in pcost calculation

diff --git a/src/glpios09.c b/src/glpios09.c
index 727761c..255e681 100644
--- a/src/glpios09.c
+++ b/src/glpios09.c
@@ -447,12 +447,19 @@ static double eval_degrad(glp_tree *T, int j, double bnd)
       double degrad;
       /* prepare lp */
       if (lp == NULL)
-      {  /* the current basis must be optimal */
+      {  glp_bfcp bfcp;
+         /* the current basis must be optimal */
          xassert(glp_get_status(P) == GLP_OPT);
          /* create a copy of mip */
          lp = glp_create_prob();
          glp_copy_prob(lp, P, 0);
          csa->work = lp;
+         glp_get_bfcp(lp, &bfcp);
+         if (bfcp.type & GLP_BF_FT)
+         {  bfcp.type &= ~GLP_BF_FT;
+            bfcp.type |= GLP_BF_BG;
+         }
+         glp_set_bfcp(lp, &bfcp);
          /* factorize and store */
          glp_factorize(lp);
          csa->bfd = bfd_create_it();
schur2.patch (text/x-patch, 4.8 KB)
commit d9944e1c142a09e20ee0cf58fa660b3d09206475
Author: Chris Matrakidis <[email protected]>
Date:   Sat Jan 14 19:36:43 2017 +0200

    use bfd_reset in pcost initialisation

diff --git a/src/bfd.c b/src/bfd.c
index 89d547a..a085298 100644
--- a/src/bfd.c
+++ b/src/bfd.c
@@ -514,6 +514,16 @@ int bfd_update(BFD *bfd, int j, int len, const int ind[], const double
       return ret;
 }
 
+void bfd_reset(BFD *bfd)
+{     /* reset LP basis factorization to initial state */
+      int ret;
+      xassert(bfd->valid);
+      xassert(bfd->type == 2);
+      scfint_reset(bfd->u.scfi);
+      bfd->upd_cnt = 0;
+      return;
+}
+
 int bfd_get_count(BFD *bfd)
 {     /* determine factorization update count */
       return bfd->upd_cnt;
diff --git a/src/bfd.h b/src/bfd.h
index 3f63099..d652c55 100644
--- a/src/bfd.h
+++ b/src/bfd.h
@@ -94,6 +94,10 @@ int bfd_update(BFD *bfd, int j, int len, const int ind[], const double
       val[]);
 /* update LP basis factorization */
 
+#define bfd_reset _glp_bfd_reset
+void bfd_reset(BFD *bfd);
+/* reset LP basis factorization to initial state */
+
 #define bfd_get_count _glp_bfd_get_count
 int bfd_get_count(BFD *bfd);
 /* determine factorization update count */
diff --git a/src/bflib/scfint.c b/src/bflib/scfint.c
index 25aed12..6e572b0 100644
--- a/src/bflib/scfint.c
+++ b/src/bflib/scfint.c
@@ -135,6 +135,26 @@ int scfint_factorize(SCFINT *fi, int n, int (*col)(void *info, int j,
       return ret;
 }
 
+void scfint_reset(SCFINT *fi)
+{     /* reset SC-factorization to initial state */
+      int k;
+      int n = fi->scf.n0;
+      xassert(n > 0);
+      /* initialize SC-factorization */
+      fi->scf.n = n;
+      fi->scf.nn = 0;
+      fi->scf.rr_ref = sva_alloc_vecs(fi->scf.sva, fi->scf.ifu.n_max);
+      fi->scf.ss_ref = sva_alloc_vecs(fi->scf.sva, fi->scf.ifu.n_max);
+      fi->scf.ifu.n = 0;
+      for (k = 1; k <= n; k++)
+      {  fi->scf.pp_ind[k] = k;
+         fi->scf.pp_inv[k] = k;
+         fi->scf.qq_ind[k] = k;
+         fi->scf.qq_inv[k] = k;
+      }
+      return;
+}
+
 int scfint_update(SCFINT *fi, int upd, int j, int len, const int ind[],
       const double val[])
 {     /* update SC-factorization after replacing j-th column of A */
diff --git a/src/bflib/scfint.h b/src/bflib/scfint.h
index 67ca65d..45142a6 100644
--- a/src/bflib/scfint.h
+++ b/src/bflib/scfint.h
@@ -63,6 +63,10 @@ int scfint_factorize(SCFINT *fi, int n, int (*col)(void *info, int j,
       int ind[], double val[]), void *info);
 /* compute SC-factorization of specified matrix A */
 
+#define scfint_reset _glp_scfint_reset
+void scfint_reset(SCFINT *fi);
+/* reset SC-factorization to initial state */
+
 #define scfint_update _glp_scfint_update
 int scfint_update(SCFINT *fi, int upd, int j, int len, const int ind[],
       const double val[]);
diff --git a/src/glpios09.c b/src/glpios09.c
index 255e681..4d3f5d6 100644
--- a/src/glpios09.c
+++ b/src/glpios09.c
@@ -406,14 +406,12 @@ struct csa
          over all up_cnt[j] subproblems */
       glp_prob *work;
       /* working problem to avoid unnecessary initializations */
-      BFD *bfd;
-      /* factorization copy to avoid identical factorizations*/
       int *head;
-      /* basis header copy to avoid identical factorizations */
+      /* basis header copy to restore factorization */
       int *r_bind;
-      /* copy of row bind values to avoid identical factorizations */
+      /* copy of row bind values to restore factorization */
       int *c_bind;
-      /* copy of column bind values to avoid identical factorizations */
+      /* copy of column bind values to restore factorization */
 };
 
 void *ios_pcost_init(glp_tree *tree)
@@ -430,7 +428,6 @@ void *ios_pcost_init(glp_tree *tree)
          csa->dn_sum[j] = csa->up_sum[j] = 0.0;
       }
       csa->work = NULL;
-      csa->bfd = NULL;
       return csa;
 }
 
@@ -462,8 +459,6 @@ static double eval_degrad(glp_tree *T, int j, double bnd)
          glp_set_bfcp(lp, &bfcp);
          /* factorize and store */
          glp_factorize(lp);
-         csa->bfd = bfd_create_it();
-         bfd_copy(csa->bfd, lp->bfd);
          /* copy factorization */
          csa->head = xcalloc(1+P->m, sizeof(int));
          csa->c_bind = xcalloc(1+P->n, sizeof(int));
@@ -490,7 +485,7 @@ static double eval_degrad(glp_tree *T, int j, double bnd)
             lp->col[jjj]->dual = P->col[jjj]->dual;
             lp->col[jjj]->bind = csa->c_bind[jjj];
          }
-         bfd_copy(lp->bfd, csa->bfd);
+         bfd_reset(lp->bfd);
          lp->valid = 1;
       }
       /* fix column x[j] at specified value */
@@ -723,7 +718,6 @@ done: *_next = sel;
          glp_delete_prob(csa->work);
          csa->work = NULL;
          /* clean up stored factorization for next round */
-         bfd_delete_it(csa->bfd);
          xfree(csa->c_bind);
          xfree(csa->r_bind);
          xfree(csa->head);
lmpx.com only provides a reader for public news (NNTP) servers. It is not affiliated with the servers or forums shown here and is not responsible for the content of articles, which is written by their respective authors.