--- /dev/null
+@@ This patch is taken from https://github.com/timothy-king/glpk-cut-log and
+@@ has the following license:
+@@
+@@ GLPK is free software: you can redistribute it and/or modify it under the
+@@ terms of the GNU General Public License as published by the Free Software
+@@ Foundation, either version 3 of the License, or (at your option) any later
+@@ version.
+@@
+From cf84e3855d8c5676557daaca434b42b2fc17dc29 Mon Sep 17 00:00:00 2001
+From: Tim King <taking@cs.nyu.edu>
+Date: Sat, 14 Dec 2013 16:03:35 -0500
+Subject: [PATCH] Adding a function for returning the iteration count.
+
+Adding GLP_ICUTADDED callback:
+- Adds a new callback location that is called after every cut is added to the pool.
+- During this callback users can request a copy of this cut and query what kind of cut it is.
+- GMI cuts also support returning what row in the tableau it came from. Users can use this to replay how the cut was generated.
+
+Logging MIR Cuts
+- Changed the IOSAUX to be defined as a linear sum of input rows.
+ This is visible from the interface
+- Logging the multipliers used to create the aggregate row in MIR.
+- Logging the delta selected by final successful cmir_sep call.
+- Logging the cset selected by final successful cmir_sep call.
+- The delta and cset are not yet user visible.
+
+The extra mir logging information is now user visible.
+
+Adding GLP_ICUTSELECT callback:
+- ios_process_cuts marks the cuts as selected after turning the cut into a row.
+- This callback happens after ios_process_cuts marks and before the pool is
+ cleared by ios_clear_pool.
+
+Cleaning up the git directory to not track autogenerated files.
+
+Branch logging
+- Adds to the integer interface a callback function for when branches are made.
+ This is called after branch_on.
+- Added a public function glp_ios_branch_log.
+ This returns the structural variable for the branch as well as the value
+ branched on, and the ids of the nodes for the down and up branches.
+ If both branches are infeasible, the ids are -1 for both dn and up.
+ If one branch is infeasible (and no branch is done), this returns -1 for
+ the infeasible branch and the current node id for the up branch.
+- Each node now has a unique ordinal associated with it.
+- Added a public function to convert the id of an active node to the unique
+ node ordinal.
+- Added a callback for when a node is closed due to being linearly infeasible.
+
+Improved the mir cut logging facilities to now include:
+- the rid id used for virtual lower/upper bound selection,
+- and the subst map,
+- Changed the size of the cset array to be the same as the new arrays (m+n)
+- Fixed a bug in returning the cset selected.
+
+Fixing a bug selecting the correct branch.
+
+Added a callback for tracking row deletion.
+
+Adding notes to node freeze and revive functions.
+
+Commenting out a few debugging printfs for gomory cuts.
+
+Changing an overly aggressive assert in MIR generate when an assignment is outside of the bounds to skipping the application of the mir rule.
+
+Removing some printfs.
+
+Weakening assert to a failure.
+
+Adding a stability failure limit to simplex.
+
+Updating the installation instructions.
+
+---
+diff --git a/configure.ac b/configure.ac
+index 0141a8a..ffb9ae4 100644
+--- a/configure.ac
++++ b/configure.ac
+@@ -1,7 +1,7 @@
+ dnl Process this file with autoconf to produce a configure script
+
+ AC_INIT([GLPK], [4.52], [bug-glpk@gnu.org])
+-
++AM_INIT_AUTOMAKE([subdir-objects])
+ AC_CONFIG_SRCDIR([src/glpk.h])
+
+ AC_CONFIG_MACRO_DIR([m4])
+diff --git a/src/Makefile.am b/src/Makefile.am
+index b39efa6..2a025bc 100644
+--- a/src/Makefile.am
++++ b/src/Makefile.am
+@@ -20,7 +20,10 @@ libglpk_la_LDFLAGS = \
+ -version-info 36:0:1 \
+ -export-symbols-regex '^glp_*'
+
++
++# all of the cut log sources are listed first
+ libglpk_la_SOURCES = \
++cutlog01.c \
+ glpapi01.c \
+ glpapi02.c \
+ glpapi03.c \
+diff --git a/src/cutlog01.c b/src/cutlog01.c
+new file mode 100644
+index 0000000..9a85bb7
+--- /dev/null
++++ b/src/cutlog01.c
+@@ -0,0 +1,452 @@
++/* cutlog01.c (api extension routines) */
++
++
++#include "glpapi.h"
++#include "glpios.h"
++
++int glp_get_it_cnt(glp_prob *P){
++ if(P == NULL){
++ return 0;
++ }else{
++ return P->it_cnt;
++ }
++}
++
++int glp_ios_get_cut(glp_tree *T, int i, int* ind, double* val, int* klass, int* type, double* rhs){
++ xassert(T != NULL);
++
++ IOSCUT* cut;
++ int len;
++ IOSAIJ* aij;
++ glp_prob* prob;
++
++ if (T->reason != GLP_ICUTADDED){
++ xerror("glp_ios_get_cut: not called during cut added.\n");
++ }
++ cut = ios_find_row(T->local, i);
++ if ( cut == NULL ) {
++ xerror("glp_ios_get_cut: called with an invalid index.");
++ }
++ len = 0;
++ for(len = 0, aij = cut->ptr; aij != NULL; aij = aij->next)
++ {
++ len++;
++ if(ind != NULL){ ind[len] = aij->j; }
++ if(val != NULL){ val[len] = aij->val; }
++ }
++ if(klass != NULL){ *klass = cut->klass; }
++ if(type != NULL){ *type = cut->type; }
++ if(rhs != NULL){ *rhs = cut->rhs; }
++ return len;
++}
++
++
++IOSAUX *ios_create_aux(int n, const int rows[], const double coeffs[]){
++ IOSAUX *aux;
++ int i;
++ aux = xmalloc(sizeof(IOSAUX));
++ aux->nrows = n;
++ aux->rows = xcalloc(1+n, sizeof(int));
++ aux->mult = xcalloc(1+n, sizeof(double));
++ aux->selected = -1;
++ aux->mir_cset = NULL;
++
++ for ( i = 1; i <= n; i++)
++ { aux->rows[i] = rows[i];
++ aux->mult[i] = coeffs[i];
++ }
++
++ return aux;
++}
++
++void ios_delete_aux(IOSAUX *aux){
++ xassert(aux != NULL);
++ xfree(aux->rows);
++ xfree(aux->mult);
++ if( aux->mir_cset != NULL ){
++ xfree(aux->mir_cset);
++ xfree(aux->mir_subst);
++ xfree(aux->mir_vlb_rows);
++ xfree(aux->mir_vub_rows);
++ }
++ xfree(aux);
++ return;
++}
++
++static void cut_set_aux(IOSCUT *cut, int n,
++ const int rows[], const double coeffs[]){
++ int i;
++ xassert( cut != NULL );
++ if( cut->aux != NULL ) {
++ ios_delete_aux(cut-> aux);
++ }
++
++ cut->aux = ios_create_aux(n, rows, coeffs);
++ xassert( cut->aux->nrows == n );
++}
++
++static void cut_set_aux_mir(IOSAUX *aux, double delta, int m, int n,
++ const char cset[], const char subst[],
++ const int vlbrs[], const int vubrs[]){
++ int i;
++ xassert( aux != NULL );
++ if ( aux->mir_cset != NULL )
++ { xfree(aux->mir_cset);
++ xfree(aux->mir_subst);
++ xfree(aux->mir_vlb_rows);
++ xfree(aux->mir_vub_rows);
++ }
++
++ aux->mir_cset = xcalloc(1+n+m, sizeof(char));
++ aux->mir_subst = xcalloc(1+n+m, sizeof(char));
++ aux->mir_vlb_rows = xcalloc(1+n+m, sizeof(int));
++ aux->mir_vub_rows = xcalloc(1+n+m, sizeof(int));
++
++ for ( i = 1; i <= n+m; i++)
++ { aux->mir_cset[i] = cset[i];
++ aux->mir_subst[i] = subst[i];
++ aux->mir_vlb_rows[i] = vlbrs[i];
++ aux->mir_vub_rows[i] = vubrs[i];
++ }
++
++ aux->mir_delta = delta;
++}
++
++void ios_cut_set_aux_mir(glp_tree *T, int ord, double delta,
++ const char cset[], const char subst[],
++ const int vlbrs[], const int vubrs[]){
++ int m, n;
++ IOSCUT *cut;
++ glp_prob *mip;
++ mip = T->mip;
++ m = mip->m;
++ n = mip->n;
++ cut = ios_find_row(T->local, ord);
++ xassert(cut != NULL);
++ cut_set_aux_mir(cut->aux, delta, m, n, cset, subst, vlbrs, vubrs);
++}
++
++void ios_cut_set_single_aux(glp_tree *T, int ord, int j){
++ IOSCUT *cut;
++ cut = ios_find_row(T->local, ord);
++ xassert(cut != NULL);
++
++ /* set up arrays */
++ int ind[1+1];
++ double coeffs[1+1];
++ ind[1] = j;
++ coeffs[1] = +1.0;
++
++ /* call general procedure */
++ cut_set_aux(cut, 1, ind, coeffs);
++}
++
++void ios_cut_set_aux(glp_tree *T, int ord, int n,
++ const int rows[], const double coeffs[]){
++ IOSCUT *cut;
++ cut = ios_find_row(T->local, ord);
++ xassert(cut != NULL);
++ cut_set_aux(cut, n, rows, coeffs);
++}
++
++int glp_ios_cut_get_aux_nrows(glp_tree *tree, int ord){
++ IOSCUT *cut;
++ IOSAUX *aux;
++ if (tree->reason != GLP_ICUTADDED){
++ xerror("glp_ios_cut_get_aux_nrows: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_nrows: not called on a valid cut.\n");
++ }
++ aux = cut->aux;
++ return (aux == NULL) ? 0 : aux->nrows;
++}
++
++void glp_ios_cut_get_aux_rows(glp_tree *tree, int ord,
++ int rows[], double coeffs[]){
++ IOSCUT *cut;
++ IOSAUX *aux;
++ int j, nrows;
++ if (tree->reason != GLP_ICUTADDED){
++ xerror("glp_ios_cut_get_aux_rows: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_rows: not called on a valid cut.\n");
++ }
++ aux = cut->aux;
++ if( aux != NULL ){
++ nrows = aux->nrows;
++ for ( j = 1; j <= nrows; j++ )
++ { if ( rows != NULL ) { rows[j] = aux->rows[j]; }
++ if ( coeffs != NULL ) { coeffs[j] = aux->mult[j]; }
++ }
++ }
++ return;
++}
++
++
++void glp_ios_cut_get_mir_subst(glp_tree *tree, int ord, char subst[]);
++/* gets mir cut substition information. */
++void glp_ios_cut_get_mir_virtual_rows(glp_tree *tree, int ord,
++ int vlb[], int vub[]);
++/* gets mir cut virtual bounds rows. */
++
++void glp_ios_cut_get_mir_cset(glp_tree *tree, int ord, char *cset){
++ glp_prob *mip;
++ IOSCUT *cut;
++ IOSAUX *aux;
++ int j, n, m;
++ if ( tree->reason != GLP_ICUTADDED ){
++ xerror("glp_ios_cut_get_aux_mir: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a cut.\n");
++ }
++ if ( cut->klass != GLP_RF_MIR ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a mir cut.\n");
++ }
++ aux = cut->aux;
++ mip = tree->mip;
++ m = mip->m;
++ n = mip->n;
++
++ if( cset != NULL ){
++ for ( j=1; j <= n+m; j++ ){
++ cset[j] = (aux->mir_cset == NULL) ? 0 : aux->mir_cset[j];
++ }
++ }
++}
++void glp_ios_cut_get_mir_subst(glp_tree *tree, int ord, char *subst){
++ glp_prob *mip;
++ IOSCUT *cut;
++ IOSAUX *aux;
++ int j, n, m;
++ if ( tree->reason != GLP_ICUTADDED ){
++ xerror("glp_ios_cut_get_aux_mir: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a cut.\n");
++ }
++ if ( cut->klass != GLP_RF_MIR ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a mir cut.\n");
++ }
++ aux = cut->aux;
++ mip = tree->mip;
++ m = mip->m;
++ n = mip->n;
++
++ if( subst != NULL ){
++ for ( j=1; j <= n+m; j++ ){
++ subst[j] = (aux->mir_subst == NULL) ? 0 : aux->mir_subst[j];
++ }
++ }
++}
++void glp_ios_cut_get_mir_virtual_rows(glp_tree *tree, int ord, int vlb_rows[], int vub_rows[]){
++ glp_prob *mip;
++ IOSCUT *cut;
++ IOSAUX *aux;
++ int j, n, m;
++ if ( tree->reason != GLP_ICUTADDED ){
++ xerror("glp_ios_cut_get_aux_mir: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a cut.\n");
++ }
++ if ( cut->klass != GLP_RF_MIR ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a mir cut.\n");
++ }
++ aux = cut->aux;
++ mip = tree->mip;
++ m = mip->m;
++ n = mip->n;
++
++ for ( j=1; j <= n+m; j++ ){
++ vlb_rows[j] = (aux->mir_vlb_rows == NULL) ? 0 : aux->mir_vlb_rows[j];
++ vub_rows[j] = (aux->mir_vub_rows == NULL) ? 0 : aux->mir_vub_rows[j];
++ }
++}
++double glp_ios_cut_get_mir_delta(glp_tree *tree, int ord){
++ glp_prob *mip;
++ IOSCUT *cut;
++ IOSAUX *aux;
++ int j, n, m;
++ if ( tree->reason != GLP_ICUTADDED ){
++ xerror("glp_ios_cut_get_aux_mir: not called during cut added.\n");
++ }
++ cut = ios_find_row(tree->local, ord);
++ if ( cut == NULL ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a cut.\n");
++ }
++ if ( cut->klass != GLP_RF_MIR ){
++ xerror("glp_ios_cut_get_aux_mir: not called on a mir cut.\n");
++ }
++ aux = cut->aux;
++ return aux->mir_delta;
++}
++
++
++void ios_cut_set_selected(IOSCUT *cut, int sel){
++#ifdef CUT_DEBUG
++ static int i = 0;
++ ++i;
++ printf("ios_cut_set_selected: %d %d %p\n", i, sel, cut);
++#endif
++
++ IOSAUX *aux;
++ aux = cut->aux;
++ if ( aux != NULL ){
++ aux->selected = sel;
++ }
++}
++
++int glp_ios_selected_cuts(glp_tree *tree, int ords[], int sel[]){
++ int len, j, N, s;
++ IOSPOOL* pool;
++ IOSCUT* cut;
++ IOSAUX* aux;
++ if ( tree == NULL ){
++ xerror("glp_ios_selected_cuts: not called with a valid tree.\n");
++ }
++ if ( tree->reason != GLP_ICUTSELECT ){
++ xerror("glp_ios_selected_cuts: not called during cut selected.\n");
++ }
++
++ pool = tree->local;
++ if ( pool == NULL ){
++ xerror("glp_ios_selected_cuts: called on a malformed tree.\n");
++ }
++
++ for (len = 0, j = 1, cut = pool->head; cut != NULL; cut = cut->next, j++)
++ { aux = cut->aux;
++#ifdef CUT_DEBUG
++ printf("glp_ios_selected_cuts: %d %p\n", j, cut);
++#endif
++ if ( aux != NULL )
++ { s = aux->selected;
++ if ( s >= 0 )
++ { len++;
++ if (ords != NULL) { ords[len] = j; }
++ if (sel != NULL) { sel[len] = s; }
++ }
++ }
++ }
++ return len;
++}
++
++int glp_ios_branch_log(glp_tree *tree, double *br_val, int* parent, int* dn, int* up){
++ IOSNPD *node;
++ int br_result, br_var;
++ int p, d, u;
++ double v;
++ glp_prob *mip;
++ if ( tree == NULL ){
++ xerror("glp_ios_branch_log: not called with a valid tree \n");
++ }
++ if ( tree->reason != GLP_LI_BRANCH ){
++ xerror("glp_ios_branch_log: not called during GLP_LI_BRANCH \n");
++ }
++ mip = tree->mip;
++ if ( mip == NULL ){
++ xerror("glp_ios_branch_log: not called with a valid tree\n");
++ }
++ br_result = tree->br_result;
++ br_var = tree->br_var;
++ switch(br_result){
++ case 0:
++ p = tree->br_node;
++ node = tree->slot[p].node;
++ break;
++ case 1:
++ case 2:
++ node = tree->curr;
++ p = node->p;
++ break;
++ default:
++ xerror("glp_ios_branch_log: br_result is not properly set \n");
++ }
++ if( node == NULL ){
++ xerror("glp_ios_branch_log: not called with a valid tree \n");
++ }
++ switch(br_result){
++ case 0:
++ v = node->br_val;
++ d = tree->dn_child;
++ u = tree->up_child;
++ break;
++ case 1:
++ v = mip->col[br_var]->prim;
++ if(tree->br_to_up){
++ d = -1;
++ u = p;
++ }else{
++ d = p;
++ u = -1;
++ }
++ break;
++ case 2:
++ v = mip->col[br_var]->prim;
++ d = -1;
++ u = -1;
++ break;
++ default:
++ xerror("glp_ios_branch_log: not called with a valid tree \n");
++ }
++
++ if(br_val != NULL){ *br_val = v; }
++ if(parent != NULL){ *parent = p; }
++ if(dn != NULL){ *dn = d; }
++ if(up != NULL){ *up = u; }
++
++ return br_var;
++}
++
++int glp_ios_node_ord(glp_tree *tree, int p){
++ IOSNPD *node;
++ if ( tree == NULL ){
++ xerror("glp_ios_node_ord: not called with a valid tree.\n");
++ }
++ if (!(1 <= p && p <= tree->nslots)){
++ xerror("glp_ios_node_ord: not called with a valid p.\n");
++ }
++ node = tree->slot[p].node;
++ return node->ord;
++}
++
++void ios_cb_rows_deleted(glp_tree *T, int n, const int* rows){
++ if (T->parm->cb_func != NULL)
++ {
++ xassert(T->reason == 0);
++ xassert(T->deleting_rows == NULL);
++ xassert(T->num_deleting_rows == 0);
++ T->num_deleting_rows = n;
++ T->deleting_rows = rows;
++
++ T->reason = GLP_LI_DELROW;
++ T->parm->cb_func(T, T->parm->cb_info);
++ T->reason = 0;
++ T->num_deleting_rows = 0;
++ T->deleting_rows = NULL;
++ }
++}
++
++int glp_ios_rows_deleted(glp_tree *tree, int* rows){
++ if ( tree == NULL ){
++ xerror("glp_ios_rows_deleted: not called with a valid tree.\n");
++ }
++ if ( tree->reason != GLP_LI_DELROW ){
++ xerror("glp_ios_rows_deleted: not called with a valid reason.\n");
++ }
++
++ int j;
++ if(rows != NULL){
++ for(j=1; j <= tree->num_deleting_rows; j++){
++ rows[j] = tree->deleting_rows[j];
++ }
++ }
++ return tree->num_deleting_rows;
++}
+diff --git a/src/glpapi06.c b/src/glpapi06.c
+index 743d6be..05d0498 100644
+--- a/src/glpapi06.c
++++ b/src/glpapi06.c
+@@ -488,6 +488,7 @@ void glp_init_smcp(glp_smcp *parm)
+ parm->out_frq = 500;
+ parm->out_dly = 0;
+ parm->presolve = GLP_OFF;
++ parm->stability_lmt = 200;
+ return;
+ }
+
+diff --git a/src/glpapi13.c b/src/glpapi13.c
+index 926dda4..55adf44 100644
+--- a/src/glpapi13.c
++++ b/src/glpapi13.c
+@@ -453,8 +453,14 @@ void glp_ios_row_attr(glp_tree *tree, int i, glp_attr *attr)
+
+ int glp_ios_pool_size(glp_tree *tree)
+ { /* determine current size of the cut pool */
+- if (tree->reason != GLP_ICUTGEN)
+- xerror("glp_ios_pool_size: operation not allowed\n");
++ switch(tree->reason)
++ { case GLP_ICUTGEN:
++ case GLP_ICUTADDED:
++ case GLP_ICUTSELECT:
++ break;
++ default:
++ xerror("glp_ios_pool_size: operation not allowed\n");
++ }
+ xassert(tree->local != NULL);
+ return tree->local->size;
+ }
+diff --git a/src/glpios.h b/src/glpios.h
+index 9e2d6b2..3b31901 100644
+--- a/src/glpios.h
++++ b/src/glpios.h
+@@ -36,6 +36,9 @@ typedef struct IOSAIJ IOSAIJ;
+ typedef struct IOSPOOL IOSPOOL;
+ typedef struct IOSCUT IOSCUT;
+
++
++typedef struct IOSAUX IOSAUX;
++
+ struct glp_tree
+ { /* branch-and-bound tree */
+ int magic;
+@@ -217,6 +220,19 @@ struct glp_tree
+ GLP_NO_BRNCH - use general selection technique */
+ int child;
+ /* subproblem reference number corresponding to br_sel */
++
++ /* start of cut log extras */
++ int br_result;
++ int br_to_up;
++ int br_node;
++ /* subproblem reference number for the just branched node.*/
++ int dn_child;
++ /* subproblem reference number for the just created down node */
++ int up_child;
++ /* subproblem reference number for the just created up node */
++
++ const int* deleting_rows;
++ int num_deleting_rows;
+ };
+
+ struct IOSLOT
+@@ -229,6 +245,8 @@ struct IOSLOT
+
+ struct IOSNPD
+ { /* node subproblem descriptor */
++ int ord;
++ /* this is a unique ordinal for each subproblem */
+ int p;
+ /* subproblem reference number (it is the index to corresponding
+ slot, i.e. slot[p] points to this descriptor) */
+@@ -393,12 +411,51 @@ struct IOSCUT
+ GLP_FX: sum a[j] * x[j] = b */
+ double rhs;
+ /* cut right-hand side */
++
++ IOSAUX *aux;
++ /* cut auxillary source information */
++
+ IOSCUT *prev;
+ /* pointer to previous cut */
+ IOSCUT *next;
+ /* pointer to next cut */
+ };
+
++struct IOSAUX
++{
++ /* aux (auxillary source information for each cut)
++ * Each cut operates on a sum of rows
++ * Each cut operates on a row that can be described using
++ * a current row or a previous cut.
++ * To generalize, we assume each row is a sum of two rows:
++ * row[r]* r_mult + pool[c] * c_mult
++ */
++ int nrows;
++ int *rows;
++ double *mult;
++
++ int selected;
++ /* when < 0 this has not yet been turned into a row
++ when >=0 this is the id of the row added. */
++
++ double mir_delta;
++ /* the delta selected by a mir cut */
++
++ char *mir_cset;
++ /* complimented set */
++ /* if this is NULL, it is implicitly all 0s */
++
++ char *mir_subst;
++ /* the substition vectors */
++
++ int *mir_vlb_rows;
++ /* the substition vectors */
++
++ int *mir_vub_rows;
++ /* the substition vectors */
++};
++
++
+ #define ios_create_tree _glp_ios_create_tree
+ glp_tree *ios_create_tree(glp_prob *mip, const glp_iocp *parm);
+ /* create branch-and-bound tree */
+@@ -615,6 +672,31 @@ int ios_choose_node(glp_tree *T);
+ int ios_choose_var(glp_tree *T, int *next);
+ /* select variable to branch on */
+
++/* functions added to retrieve information */
++IOSAUX *ios_create_aux(int n, const int rows[], const double coeffs[]);
++/* creates an aux with n rows */
++
++void ios_delete_aux(IOSAUX *aux);
++/* deletes an aux */
++
++void ios_cut_set_single_aux(glp_tree *T, int ord, int j);
++/* sets the aux of row ord to be a single row j */
++
++
++void ios_cut_set_aux(glp_tree *T, int ord, int n,
++ const int rows[], const double coeffs[]);
++/* sets an arbitrary aux sum */
++
++void ios_cut_set_aux_mir(glp_tree *T, int ord, double delta,
++ const char cset[], const char subst[],
++ const int vlbrs[], const int vubrs[]);
++/* sets the extra mir information */
++
++void ios_cut_set_selected(IOSCUT *cut, int i);
++/* the cut has been added as row i */
++
++void ios_cb_rows_deleted(glp_tree *T, int n, const int* rows);
++
+ #endif
+
+ /* eof */
+diff --git a/src/glpios01.c b/src/glpios01.c
+index 70798fd..d9e5703 100644
+--- a/src/glpios01.c
++++ b/src/glpios01.c
+@@ -159,6 +159,16 @@ glp_tree *ios_create_tree(glp_prob *mip, const glp_iocp *parm)
+ tree->stop = 0;
+ /* create the root subproblem, which initially is identical to
+ the original MIP */
++
++ tree->br_result = 0;
++ tree->br_to_up = 0;
++ tree->br_node = 0;
++ tree->dn_child = 0;
++ tree->up_child = 0;
++
++ tree->deleting_rows = NULL;
++ tree->num_deleting_rows = 0;
++
+ new_node(tree, NULL);
+ return tree;
+ }
+@@ -278,6 +288,7 @@ void ios_revive_node(glp_tree *tree, int p)
+ double *val;
+ ind = xcalloc(1+n, sizeof(int));
+ val = xcalloc(1+n, sizeof(double));
++ /* maintains the row order during revival */
+ for (r = node->r_ptr; r != NULL; r = r->next)
+ { i = glp_add_rows(mip, 1);
+ glp_set_row_name(mip, i, r->name);
+@@ -457,6 +468,13 @@ void ios_freeze_node(glp_tree *tree)
+ double *val;
+ ind = xcalloc(1+n, sizeof(int));
+ val = xcalloc(1+n, sizeof(double));
++ /* Added rows are stored in the same order!
++ * This is done by 2 reversals.
++ * - Going from i = m down pred_m.
++ * - Storing the list by a stack "push"
++ * node->r_ptr = r;
++ * r->next = node->r_ptr;
++ */
+ for (i = m; i > pred_m; i--)
+ { GLPROW *row = mip->row[i];
+ IOSROW *r;
+@@ -501,6 +519,8 @@ void ios_freeze_node(glp_tree *tree)
+ xassert(nrs > 0);
+ num = xcalloc(1+nrs, sizeof(int));
+ for (i = 1; i <= nrs; i++) num[i] = root_m + i;
++ /* To not call ios_cb_rows_deleted here.
++ These rows have been saved earlier. */
+ glp_del_rows(mip, nrs, num);
+ xfree(num);
+ }
+@@ -636,6 +656,9 @@ static IOSNPD *new_node(glp_tree *tree, IOSNPD *parent)
+ tree->a_cnt++;
+ tree->n_cnt++;
+ tree->t_cnt++;
++
++ node->ord = tree->t_cnt;
++
+ /* increase the number of child subproblems */
+ if (parent == NULL)
+ xassert(p == 1);
+@@ -818,6 +841,8 @@ void ios_delete_tree(glp_tree *tree)
+ xassert(nrs > 0);
+ num = xcalloc(1+nrs, sizeof(int));
+ for (i = 1; i <= nrs; i++) num[i] = tree->orig_m + i;
++ /* Do not call ios_cb_rows_deleted here.
++ This does not help log information. */
+ glp_del_rows(mip, nrs, num);
+ xfree(num);
+ }
+@@ -1417,6 +1442,7 @@ int ios_add_row(glp_tree *tree, IOSPOOL *pool,
+ cut->rhs = rhs;
+ cut->prev = pool->tail;
+ cut->next = NULL;
++ cut->aux = NULL;
+ if (cut->prev == NULL)
+ pool->head = cut;
+ else
+@@ -1517,6 +1543,11 @@ void ios_del_row(glp_tree *tree, IOSPOOL *pool, int i)
+ cut->ptr = aij->next;
+ dmp_free_atom(tree->pool, aij, sizeof(IOSAIJ));
+ }
++
++ if (cut->aux != NULL){
++ ios_delete_aux(cut->aux);
++ }
++
+ dmp_free_atom(tree->pool, cut, sizeof(IOSCUT));
+ pool->size--;
+ return;
+@@ -1535,6 +1566,10 @@ void ios_clear_pool(glp_tree *tree, IOSPOOL *pool)
+ cut->ptr = aij->next;
+ dmp_free_atom(tree->pool, aij, sizeof(IOSAIJ));
+ }
++
++ if (cut->aux != NULL){
++ ios_delete_aux(cut->aux);
++ }
+ dmp_free_atom(tree->pool, cut, sizeof(IOSCUT));
+ }
+ pool->size = 0;
+diff --git a/src/glpios03.c b/src/glpios03.c
+index 80d701b..03e9208 100644
+--- a/src/glpios03.c
++++ b/src/glpios03.c
+@@ -358,12 +358,12 @@ static void fix_by_red_cost(glp_tree *T)
+ * 2 - both branches are hopeless and have been pruned; new subproblem
+ * selection is needed to continue the search. */
+
+-static int branch_on(glp_tree *T, int j, int next)
++static int branch_on(glp_tree *T, int j, int next, int clone[], int* to_up)
+ { glp_prob *mip = T->mip;
+ IOSNPD *node;
+ int m = mip->m;
+ int n = mip->n;
+- int type, dn_type, up_type, dn_bad, up_bad, p, ret, clone[1+2];
++ int type, dn_type, up_type, dn_bad, up_bad, p, ret;
+ double lb, ub, beta, new_ub, new_lb, dn_lp, up_lp, dn_bnd, up_bnd;
+ /* determine bounds and value of x[j] in optimal solution to LP
+ relaxation of the current subproblem */
+@@ -431,6 +431,7 @@ static int branch_on(glp_tree *T, int j, int next)
+ else
+ xassert(mip != mip);
+ ret = 1;
++ *to_up = 0; /* up is bad. Do not go to up. */
+ goto done;
+ }
+ else if (dn_bad)
+@@ -449,6 +450,7 @@ static int branch_on(glp_tree *T, int j, int next)
+ else
+ xassert(mip != mip);
+ ret = 1;
++ *to_up = 1; /* down is bad. Go to up. */
+ goto done;
+ }
+ /* both down- and up-branches seem to be hopeful */
+@@ -702,7 +704,8 @@ static void remove_cuts(glp_tree *T)
+ }
+ }
+ if (cnt > 0)
+- { glp_del_rows(T->mip, cnt, num);
++ { ios_cb_rows_deleted(T, cnt, num);
++ glp_del_rows(T->mip, cnt, num);
+ #if 0
+ xprintf("%d inactive cut(s) removed\n", cnt);
+ #endif
+@@ -784,6 +787,7 @@ static void display_cut_info(glp_tree *T)
+
+ int ios_driver(glp_tree *T)
+ { int p, curr_p, p_stat, d_stat, ret;
++ int branch_clones[1 + 2];
+ #if 1 /* carry out to glp_tree */
+ int pred_p = 0;
+ /* if the current subproblem has been just created due to
+@@ -1010,6 +1014,12 @@ more: /* minor loop starts here */
+ if (T->parm->msg_lev >= GLP_MSG_DBG)
+ xprintf("LP relaxation has no feasible solution\n");
+ /* prune the branch */
++ if (T->parm->cb_func != NULL)
++ { xassert(T->reason == 0);
++ T->reason = GLP_LI_CLOSE;
++ T->parm->cb_func(T, T->parm->cb_info);
++ T->reason = 0;
++ }
+ goto fath;
+ }
+ else
+@@ -1219,6 +1229,14 @@ more: /* minor loop starts here */
+ ios_process_cuts(T);
+ T->reason = 0;
+ }
++ /* if the local cut pool is not empty and the callback func is there,
++ this gives the callback the chance to see what was selected. */
++ if (T->parm->cb_func != NULL && T->local->size > 0)
++ { xassert(T->reason == 0);
++ T->reason = GLP_ICUTSELECT;
++ T->parm->cb_func(T, T->parm->cb_info);
++ T->reason = 0;
++ }
+ /* clear the local cut pool */
+ ios_clear_pool(T, T->local);
+ /* perform re-optimization, if necessary */
+@@ -1255,7 +1273,34 @@ more: /* minor loop starts here */
+ T->br_var = ios_choose_var(T, &T->br_sel);
+ /* perform actual branching */
+ curr_p = T->curr->p;
+- ret = branch_on(T, T->br_var, T->br_sel);
++ xassert(T->br_to_up == 0);
++ ret = branch_on(T, T->br_var, T->br_sel, branch_clones, &T->br_to_up);
++ if (T->parm->cb_func != NULL)
++ { xassert(T->reason == 0);
++ xassert(T->br_node == 0);
++ xassert(T->dn_child == 0);
++ xassert(T->up_child == 0);
++ // record a branch here
++ T->reason = GLP_LI_BRANCH;
++ // at this point T->br_var is the branching variable
++ T->br_node = curr_p;
++ T->br_result = ret;
++ if(ret == 0){
++ T->dn_child = branch_clones[1];
++ T->up_child = branch_clones[2];
++ }
++ T->parm->cb_func(T, T->parm->cb_info);
++ T->reason = 0;
++ T->br_node = 0;
++ T->dn_child = 0;
++ T->up_child = 0;
++ if (T->stop)
++ { ret = GLP_ESTOP;
++ goto done;
++ }
++ }
++ T->br_to_up = 0;
++
+ T->br_var = T->br_sel = 0;
+ if (ret == 0)
+ { /* both branches have been created */
+diff --git a/src/glpios05.c b/src/glpios05.c
+index b9322b9..caf69d7 100644
+--- a/src/glpios05.c
++++ b/src/glpios05.c
+@@ -65,6 +65,9 @@ static void gen_cut(glp_tree *tree, struct worka *worka, int j)
+ double *phi = worka->phi;
+ int i, k, len, kind, stat;
+ double lb, ub, alfa, beta, ksi, phi1, rhs;
++ int input_j;
++ input_j = j;
++
+ /* compute row of the simplex tableau, which (row) corresponds
+ to specified basic variable xB[i] = x[m+j]; see (23) */
+ len = glp_eval_tab_row(mip, m+j, ind, val);
+@@ -73,6 +76,7 @@ static void gen_cut(glp_tree *tree, struct worka *worka, int j)
+ if it would be computed with formula (27); it is assumed that
+ beta[i] is fractional enough */
+ beta = mip->col[j]->prim;
++
+ /* compute cut coefficients phi and right-hand side rho, which
+ correspond to formula (30); dense format is used, because rows
+ of the simplex tableau is usually dense */
+@@ -104,6 +108,7 @@ static void gen_cut(glp_tree *tree, struct worka *worka, int j)
+ xassert(stat != GLP_BS);
+ /* determine row coefficient ksi[i,j] at xN[j]; see (23) */
+ ksi = val[j];
++ /* printf("%d %4.15f ", k, ksi); */
+ /* if ksi[i,j] is too large in the magnitude, do not generate
+ the cut */
+ if (fabs(ksi) > 1e+05) goto fini;
+@@ -135,6 +140,7 @@ static void gen_cut(glp_tree *tree, struct worka *worka, int j)
+ /* y[j] is integer */
+ if (fabs(alfa - floor(alfa + 0.5)) < 1e-10)
+ { /* alfa[i,j] is close to nearest integer; skip it */
++ /* printf("(skip)"); */
+ goto skip;
+ }
+ else if (f(alfa) <= f(beta))
+@@ -170,6 +176,13 @@ static void gen_cut(glp_tree *tree, struct worka *worka, int j)
+ }
+ skip: ;
+ }
++ /* printf("\n"); */
++ /* for (i = 1; i <= m+n; i++) */
++ /* { */
++ /* printf("%i %f, ", i, phi[i]); */
++ /* } */
++ /* printf("\n"); */
++
+ /* now the cut has the form sum_k phi[k] * x[k] >= rho, where cut
+ coefficients are stored in the array phi in dense format;
+ x[1,...,m] are auxiliary variables, x[m+1,...,m+n] are struc-
+@@ -218,8 +231,20 @@ skip: ;
+ ios_add_cut_row(tree, pool, GLP_RF_GMI, len, ind, val, GLP_LO,
+ rhs);
+ #else
+- glp_ios_add_row(tree, NULL, GLP_RF_GMI, 0, len, ind, val, GLP_LO,
+- rhs);
++ int ord;
++ ord = glp_ios_add_row(tree, NULL, GLP_RF_GMI, 0, len, ind, val,
++ GLP_LO, rhs);
++ ios_cut_set_single_aux(tree, ord, input_j);
++
++ /* printf("ord: % d beta %f\n", ord, beta); */
++
++ /** callback for a cut being added to the cut pool */
++ if (tree->parm->cb_func != NULL)
++ { xassert(tree->reason == GLP_ICUTGEN);
++ tree->reason = GLP_ICUTADDED;
++ tree->parm->cb_func(tree, tree->parm->cb_info);
++ tree->reason = GLP_ICUTGEN;
++ }
+ #endif
+ fini: return;
+ }
+diff --git a/src/glpios06.c b/src/glpios06.c
+index 2bd0453..4bd3cf5 100644
+--- a/src/glpios06.c
++++ b/src/glpios06.c
+@@ -100,6 +100,29 @@ struct MIR
+ /* sparse vector of cutting plane coefficients, alpha[k] */
+ double cut_rhs;
+ /* right-hand size of the cutting plane, beta */
++
++ /*-------------------------------------------------------------*/
++ /* Extras I've added to reproduce a cut externally */
++ double cut_delta;
++ /* the delta used for generating the cut */
++ char *cut_cset; /* char cut_vec[1+m+n]; */
++ /* cut_cset[k], 1 <= k <= m+n, is set to true if structural
++ variable x[k] was complemented in the cut:
++ 0 - x[k] has been not been complemented
++ non 0 - x[k] has been complemented */
++ int *vlb_rows; /* int vlb_rows[1+m+n]; */
++ /* vlb_rows[k], 1 <= k <= m+n,
++ * vlb_rows[k] <= 0 if virtual lower bound has not been set
++ * vlb_rows[k] = r if virtual lower bound was set using the row r
++ */
++ int *vub_rows; /* int vub_rows[1+m+n]; */
++ /* vub_rows[k], 1 <= k <= m+n,
++ * vub_rows[k] <= 0 if virtual upper bound has not been set
++ * vub_rows[k] = r if virtual upper bound was set using the row r
++ */
++
++ double *agg_coeffs;
++ /* coefficients used to multiply agg_coeffs */
+ };
+
+ /***********************************************************************
+@@ -146,6 +169,7 @@ static void set_row_attrib(glp_tree *tree, struct MIR *mir)
+ xassert(row != row);
+ }
+ mir->vlb[k] = mir->vub[k] = 0;
++ mir->vlb_rows[k] = mir->vub_rows[k] = 0;
+ }
+ return;
+ }
+@@ -181,6 +205,7 @@ static void set_col_attrib(glp_tree *tree, struct MIR *mir)
+ xassert(col != col);
+ }
+ mir->vlb[k] = mir->vub[k] = 0;
++ mir->vlb_rows[k] = mir->vub_rows[k] = 0;
+ }
+ return;
+ }
+@@ -230,6 +255,7 @@ static void set_var_bounds(glp_tree *tree, struct MIR *mir)
+ { /* set variable lower bound for x1 */
+ mir->lb[k1] = - a2 / a1;
+ mir->vlb[k1] = k2;
++ mir->vlb_rows[k1] = i;
+ /* the row should not be used */
+ mir->skip[i] = 1;
+ }
+@@ -240,6 +266,7 @@ static void set_var_bounds(glp_tree *tree, struct MIR *mir)
+ { /* set variable upper bound for x1 */
+ mir->ub[k1] = - a2 / a1;
+ mir->vub[k1] = k2;
++ mir->vub_rows[k1] = i;
+ /* the row should not be used */
+ mir->skip[i] = 1;
+ }
+@@ -313,6 +340,13 @@ void *ios_mir_init(glp_tree *tree)
+ mir->subst = xcalloc(1+m+n, sizeof(char));
+ mir->mod_vec = ios_create_vec(m+n);
+ mir->cut_vec = ios_create_vec(m+n);
++
++ /* added */
++ mir->cut_cset = xcalloc(1+m+n, sizeof(char));
++ mir->vlb_rows = xcalloc(1+m+n, sizeof(int));
++ mir->vub_rows = xcalloc(1+m+n, sizeof(int));
++ mir->agg_coeffs = xcalloc(1+MAXAGGR, sizeof(double));
++
+ /* set global row attributes */
+ set_row_attrib(tree, mir);
+ /* set global column attributes */
+@@ -405,6 +439,9 @@ static void initial_agg_row(glp_tree *tree, struct MIR *mir, int i)
+ mir->skip[i] = 2;
+ mir->agg_cnt = 1;
+ mir->agg_row[1] = i;
++
++ mir->agg_coeffs[1] = +1.0;
++
+ /* use x[i] - sum a[i,j] * x[m+j] = 0, where x[i] is auxiliary
+ variable of row i, x[m+j] are structural variables */
+ ios_clear_vec(mir->agg_vec);
+@@ -784,19 +821,19 @@ static int cmir_cmp(const void *p1, const void *p2)
+
+ static double cmir_sep(const int n, const double a[], const double b,
+ const double u[], const double x[], const double s,
+- double alpha[], double *beta, double *gamma)
++ double alpha[], double *beta, double *gamma,
++ double* delta, char cset[])
+ { int fail, j, k, nv, v;
+- double delta, eps, d_try[1+3], r, r_best;
+- char *cset;
++ double eps, d_try[1+3], r, r_best;
+ struct vset *vset;
+ /* allocate working arrays */
+- cset = xcalloc(1+n, sizeof(char));
++ //cset = xcalloc(1+n, sizeof(char));
+ vset = xcalloc(1+n, sizeof(struct vset));
+ /* choose initial C */
+ for (j = 1; j <= n; j++)
+ cset[j] = (char)(x[j] >= 0.5 * u[j]);
+ /* choose initial delta */
+- r_best = delta = 0.0;
++ r_best = (*delta) = 0.0;
+ for (j = 1; j <= n; j++)
+ { xassert(a[j] != 0.0);
+ /* if x[j] is close to its bounds, skip it */
+@@ -809,16 +846,16 @@ static double cmir_sep(const int n, const double a[], const double b,
+ /* compute violation */
+ r = - (*beta) - (*gamma) * s;
+ for (k = 1; k <= n; k++) r += alpha[k] * x[k];
+- if (r_best < r) r_best = r, delta = fabs(a[j]);
++ if (r_best < r) r_best = r, (*delta) = fabs(a[j]);
+ }
+ if (r_best < 0.001) r_best = 0.0;
+ if (r_best == 0.0) goto done;
+- xassert(delta > 0.0);
++ xassert((*delta) > 0.0);
+ /* try to increase violation by dividing delta by 2, 4, and 8,
+ respectively */
+- d_try[1] = delta / 2.0;
+- d_try[2] = delta / 4.0;
+- d_try[3] = delta / 8.0;
++ d_try[1] = (*delta) / 2.0;
++ d_try[2] = (*delta) / 4.0;
++ d_try[3] = (*delta) / 8.0;
+ for (j = 1; j <= 3; j++)
+ { /* construct c-MIR inequality */
+ fail = cmir_ineq(n, a, b, u, cset, d_try[j], alpha, beta,
+@@ -827,7 +864,7 @@ static double cmir_sep(const int n, const double a[], const double b,
+ /* compute violation */
+ r = - (*beta) - (*gamma) * s;
+ for (k = 1; k <= n; k++) r += alpha[k] * x[k];
+- if (r_best < r) r_best = r, delta = d_try[j];
++ if (r_best < r) r_best = r, (*delta) = d_try[j];
+ }
+ /* build subset of variables lying strictly between their bounds
+ and order it by nondecreasing values of |x[j] - u[j]/2| */
+@@ -849,7 +886,7 @@ static double cmir_sep(const int n, const double a[], const double b,
+ /* replace x[j] by its complement or vice versa */
+ cset[j] = (char)!cset[j];
+ /* construct c-MIR inequality */
+- fail = cmir_ineq(n, a, b, u, cset, delta, alpha, beta, gamma);
++ fail = cmir_ineq(n, a, b, u, cset, (*delta), alpha, beta, gamma);
+ /* restore the variable */
+ cset[j] = (char)!cset[j];
+ /* do not replace the variable in case of failure */
+@@ -860,10 +897,11 @@ static double cmir_sep(const int n, const double a[], const double b,
+ if (r_best < r) r_best = r, cset[j] = (char)!cset[j];
+ }
+ /* construct the best c-MIR inequality chosen */
+- fail = cmir_ineq(n, a, b, u, cset, delta, alpha, beta, gamma);
++ fail = cmir_ineq(n, a, b, u, cset, (*delta), alpha, beta, gamma);
+ xassert(!fail);
++
+ done: /* free working arrays */
+- xfree(cset);
++
+ xfree(vset);
+ /* return to the calling routine */
+ return r_best;
+@@ -874,7 +912,8 @@ static double generate(struct MIR *mir)
+ int m = mir->m;
+ int n = mir->n;
+ int j, k, kk, nint;
+- double s, *u, *x, *alpha, r_best = 0.0, b, beta, gamma;
++ double s, *u, *x, *alpha, r_best = 0.0, b, beta, gamma, delta;
++ char *cset;
+ ios_copy_vec(mir->cut_vec, mir->mod_vec);
+ mir->cut_rhs = mir->mod_rhs;
+ /* remove small terms, which can appear due to substitution of
+@@ -923,6 +962,7 @@ static double generate(struct MIR *mir)
+ u = xcalloc(1+nint, sizeof(double));
+ x = xcalloc(1+nint, sizeof(double));
+ alpha = xcalloc(1+nint, sizeof(double));
++ cset = xcalloc(1+nint, sizeof(char));
+ /* determine u and x */
+ for (j = 1; j <= nint; j++)
+ { k = mir->cut_vec->ind[j];
+@@ -936,6 +976,7 @@ static double generate(struct MIR *mir)
+ x[j] = mir->ub[k] - mir->x[k];
+ else
+ xassert(k != k);
++ if(!(x[j] >= -0.001)) { goto skip; }
+ xassert(x[j] >= -0.001);
+ if (x[j] < 0.0) x[j] = 0.0;
+ }
+@@ -965,6 +1006,7 @@ static double generate(struct MIR *mir)
+ }
+ else
+ xassert(k != k);
++ if(!(x >= -0.001)) { goto skip; }
+ xassert(x >= -0.001);
+ if (x < 0.0) x = 0.0;
+ s -= mir->cut_vec->val[j] * x;
+@@ -973,7 +1015,7 @@ static double generate(struct MIR *mir)
+ /* apply heuristic to obtain most violated c-MIR inequality */
+ b = mir->cut_rhs;
+ r_best = cmir_sep(nint, mir->cut_vec->val, b, u, x, s, alpha,
+- &beta, &gamma);
++ &beta, &gamma, &delta, cset);
+ if (r_best == 0.0) goto skip;
+ xassert(r_best > 0.0);
+ /* convert to raw cut */
+@@ -988,10 +1030,22 @@ static double generate(struct MIR *mir)
+ #if _MIR_DEBUG
+ ios_check_vec(mir->cut_vec);
+ #endif
++ /* added */
++ mir->cut_delta = delta;
++ // this is not a great place for resetting the array,
++ // but it should be sufficient
++ for (j = 1; j <= n+m; j++)
++ mir->cut_cset[j] = 0;
++ for (j = 1; j <= nint; j++)
++ { k = mir->cut_vec->ind[j];
++ xassert(m <= k && k <= m+n);
++ mir->cut_cset[k] = cset[j];
++ }
+ skip: /* free working arrays */
+ xfree(u);
+ xfree(x);
+ xfree(alpha);
++ xfree(cset);
+ done: return r_best;
+ }
+
+@@ -1199,8 +1253,25 @@ static void add_cut(glp_tree *tree, struct MIR *mir)
+ ios_add_cut_row(tree, pool, GLP_RF_MIR, len, ind, val, GLP_UP,
+ mir->cut_rhs);
+ #else
+- glp_ios_add_row(tree, NULL, GLP_RF_MIR, 0, len, ind, val, GLP_UP,
++ int ord;
++ ord = glp_ios_add_row(tree, NULL, GLP_RF_MIR, 0, len, ind, val, GLP_UP,
+ mir->cut_rhs);
++ ios_cut_set_aux(tree, ord, mir->agg_cnt, mir->agg_row, mir->agg_coeffs);
++ ios_cut_set_aux_mir(tree, ord, mir->cut_delta,
++ mir->cut_cset, mir->subst,
++ mir->vlb_rows, mir->vub_rows);
++
++ /** callback for a cut being added to the cut pool */
++ /* printf("mir tree parm %p %d\n", tree->parm->cb_func, ord); */
++ /* printf(" agg_rhs %f\n", mir->agg_rhs); */
++ /* printf(" mod_rhs %f\n", mir->mod_rhs); */
++ /* printf(" cut_rhs %f\n", mir->cut_rhs); */
++ if (tree->parm->cb_func != NULL)
++ { xassert(tree->reason == GLP_ICUTGEN);
++ tree->reason = GLP_ICUTADDED;
++ tree->parm->cb_func(tree, tree->parm->cb_info);
++ tree->reason = GLP_ICUTGEN;
++ }
+ #endif
+ xfree(ind);
+ xfree(val);
+@@ -1216,6 +1287,7 @@ static int aggregate_row(glp_tree *tree, struct MIR *mir)
+ IOSVEC *v;
+ int ii, j, jj, k, kk, kappa = 0, ret = 0;
+ double d1, d2, d, d_max = 0.0;
++ double guass_coeff;
+ /* choose appropriate structural variable in the aggregated row
+ to be substituted */
+ for (j = 1; j <= mir->agg_vec->nnz; j++)
+@@ -1300,8 +1372,9 @@ static int aggregate_row(glp_tree *tree, struct MIR *mir)
+ xassert(j != 0);
+ jj = v->pos[kappa];
+ xassert(jj != 0);
+- ios_linear_comb(mir->agg_vec,
+- - mir->agg_vec->val[j] / v->val[jj], v);
++ guass_coeff = - mir->agg_vec->val[j] / v->val[jj];
++ mir->agg_coeffs[mir->agg_cnt] = guass_coeff;
++ ios_linear_comb(mir->agg_vec, guass_coeff, v);
+ ios_delete_vec(v);
+ ios_set_vj(mir->agg_vec, kappa, 0.0);
+ #if _MIR_DEBUG
+@@ -1440,6 +1513,13 @@ void ios_mir_term(void *gen)
+ xfree(mir->subst);
+ ios_delete_vec(mir->mod_vec);
+ ios_delete_vec(mir->cut_vec);
++
++ /* added */
++ xfree(mir->cut_cset);
++ xfree(mir->vlb_rows);
++ xfree(mir->vub_rows);
++ xfree(mir->agg_coeffs);
++
+ xfree(mir);
+ return;
+ }
+diff --git a/src/glpios07.c b/src/glpios07.c
+index 0892550..9f742e6 100644
+--- a/src/glpios07.c
++++ b/src/glpios07.c
+@@ -543,6 +543,14 @@ void ios_cov_gen(glp_tree *tree)
+ /* add the cut to the cut pool */
+ glp_ios_add_row(tree, NULL, GLP_RF_COV, 0, len, ind, val,
+ GLP_UP, val[0]);
++
++ /** callback for a cut being added to the cut pool */
++ if (tree->parm->cb_func != NULL)
++ { xassert(tree->reason == GLP_ICUTGEN);
++ tree->reason = GLP_ICUTADDED;
++ tree->parm->cb_func(tree, tree->parm->cb_info);
++ tree->reason = GLP_ICUTGEN;
++ }
+ }
+ /* free working arrays */
+ xfree(ind);
+diff --git a/src/glpios08.c b/src/glpios08.c
+index 2d2fcd5..3aad808 100644
+--- a/src/glpios08.c
++++ b/src/glpios08.c
+@@ -129,6 +129,15 @@ void ios_clq_gen(glp_tree *T, void *G_)
+ /* add cut inequality to local cut pool */
+ glp_ios_add_row(T, NULL, GLP_RF_CLQ, 0, len, ind, val, GLP_UP,
+ rhs);
++
++ /** callback for a cut being added to the cut pool */
++ if (T->parm->cb_func != NULL)
++ { xassert(T->reason == GLP_ICUTGEN);
++ T->reason = GLP_ICUTADDED;
++ T->parm->cb_func(T, T->parm->cb_info);
++ T->reason = GLP_ICUTGEN;
++ }
++
+ skip: /* free working arrays */
+ tfree(ind);
+ tfree(val);
+diff --git a/src/glpios11.c b/src/glpios11.c
+index c40e9a5..82d5698 100644
+--- a/src/glpios11.c
++++ b/src/glpios11.c
+@@ -193,6 +193,9 @@ void ios_process_cuts(glp_tree *T)
+ glp_set_mat_row(T->mip, i, len, ind, val);
+ xassert(cut->type == GLP_LO || cut->type == GLP_UP);
+ glp_set_row_bnds(T->mip, i, cut->type, cut->rhs, cut->rhs);
++
++ /* setting this as selected */
++ ios_cut_set_selected(cut, i);
+ }
+ /* free working arrays */
+ xfree(info);
+diff --git a/src/glpk.h b/src/glpk.h
+index 75c292c..e117782 100644
+--- a/src/glpk.h
++++ b/src/glpk.h
+@@ -130,6 +130,7 @@ typedef struct
+ int out_frq; /* spx.out_frq */
+ int out_dly; /* spx.out_dly (milliseconds) */
+ int presolve; /* enable/disable using LP presolver */
++ int stability_lmt; /* maximum number of check stability failures before stopping */
+ double foo_bar[36]; /* (reserved) */
+ } glp_smcp;
+
+@@ -229,6 +230,11 @@ typedef struct
+ #define GLP_IBRANCH 0x05 /* request for branching */
+ #define GLP_ISELECT 0x06 /* request for subproblem selection */
+ #define GLP_IPREPRO 0x07 /* request for preprocessing */
++#define GLP_ICUTADDED 0x08 /* cut was added to the pool */
++#define GLP_ICUTSELECT 0x09 /* cuts were selected as rows */
++#define GLP_LI_BRANCH 0x10 /* a branch was made */
++#define GLP_LI_CLOSE 0x11 /* an active node was closed */
++#define GLP_LI_DELROW 0x12 /* an active node was closed */
+
+ /* branch selection indicator: */
+ #define GLP_NO_BRNCH 0 /* select no branch */
+@@ -1057,6 +1063,54 @@ int glp_top_sort(glp_graph *G, int v_num);
+ int glp_wclique_exact(glp_graph *G, int v_wgt, double *sol, int v_set);
+ /* find maximum weight clique with exact algorithm */
+
++/*******************************************/
++/*** CUT LOG ***/
++/*******************************************/
++
++int glp_get_it_cnt(glp_prob *P);
++/* get the iteration count of the current problem */
++
++
++int glp_ios_get_cut(glp_tree *T, int i, int ind[], double val[], int* klass, int* type, double* rhs);
++/* determine reason for calling the callback routine */
++
++int glp_ios_cut_get_aux_nrows(glp_tree *tree, int ord);
++/* gets the number of rows used to generate a cut. */
++
++
++void glp_ios_cut_get_aux_rows(glp_tree *tree, int ord,
++ int rows[], double coeffs[]);
++/* gets a cut as an input sequence of rows times coefficients. */
++
++
++void glp_ios_cut_get_mir_cset(glp_tree *tree, int ord, char cset[]);
++/* gets mir cut complement set. */
++double glp_ios_cut_get_mir_delta(glp_tree *tree, int ord);
++/* gets mir cut delta. */
++void glp_ios_cut_get_mir_subst(glp_tree *tree, int ord, char subst[]);
++/* gets mir cut substition information. */
++void glp_ios_cut_get_mir_virtual_rows(glp_tree *tree, int ord,
++ int vlb_rows[], int vub_rows[]);
++/* gets mir cut virtual bounds rows. */
++
++int glp_ios_selected_cuts(glp_tree *tree, int ords[], int sel[]);
++/* gets the list of selected cuts.
++ Can only be called when GLP_ICUTSELECT */
++
++
++int glp_ios_branch_log(glp_tree *tree, double *val, int* parent, int* dn, int* up);
++/* can only be called when GLP_LI_BRANCH.
++ * If id is non-null, returns the id of the structural variable branched upon.
++ * If val is non-null, it is set to the value branched upon.
++ * If parent is non-null, it is set to node id of the node branched upon.
++ * If dn is non-null, it is set to node id of the newly created down node.
++ * If up is non-null, it is set to node id of the newly created up node.
++ */
++
++int glp_ios_node_ord(glp_tree *tree, int node_p);
++
++int glp_ios_rows_deleted(glp_tree *tree, int* rows);
++
+ #ifdef __cplusplus
+ }
+ #endif
+diff --git a/src/glpspx01.c b/src/glpspx01.c
+index 5c17114..caaab27 100644
+--- a/src/glpspx01.c
++++ b/src/glpspx01.c
+@@ -238,6 +238,9 @@ struct csa
+ double *work2; /* double work2[1+m]; */
+ double *work3; /* double work3[1+m]; */
+ double *work4; /* double work4[1+m]; */
++
++ /** Things Tim has added. */
++ int stability_failures;
+ };
+
+ static const double kappa = 0.10;
+@@ -400,6 +403,8 @@ static void init_csa(struct csa *csa, glp_prob *lp)
+ csa->refct = 0;
+ memset(&refsp[1], 0, (m+n) * sizeof(char));
+ for (j = 1; j <= n; j++) gamma[j] = 1.0;
++
++ csa->stability_failures = 0;
+ return;
+ }
+
+@@ -2647,6 +2652,11 @@ loop: /* main loop starts here */
+ if (check_stab(csa, parm->tol_bnd))
+ { /* there are excessive bound violations due to round-off
+ errors */
++ csa->stability_failures++;
++ if (csa->stability_failures >= parm->stability_lmt){
++ ret = GLP_EINSTAB;
++ goto done;
++ }
+ if (parm->msg_lev >= GLP_MSG_ERR)
+ xprintf("Warning: numerical instability (primal simplex,"
+ " phase %s)\n", csa->phase == 1 ? "I" : "II");
+diff --git a/src/glpspx02.c b/src/glpspx02.c
+index 19152a6..de3d677 100644
+--- a/src/glpspx02.c
++++ b/src/glpspx02.c
+@@ -288,6 +288,9 @@ struct csa
+ double *work2; /* double work2[1+m]; */
+ double *work3; /* double work3[1+m]; */
+ double *work4; /* double work4[1+m]; */
++
++ /* count the number of stability failures */
++ int stability_failures;
+ };
+
+ static const double kappa = 0.10;
+@@ -501,6 +504,8 @@ static void init_csa(struct csa *csa, glp_prob *lp)
+ csa->refct = 0;
+ memset(&refsp[1], 0, (m+n) * sizeof(char));
+ for (i = 1; i <= m; i++) gamma[i] = 1.0;
++
++ csa->stability_failures = 0;
+ return;
+ }
+
+@@ -2741,7 +2746,13 @@ loop: /* main loop starts here */
+ /* make sure that the current basic solution remains dual
+ feasible */
+ if (check_stab(csa, parm->tol_dj) != 0)
+- { if (parm->msg_lev >= GLP_MSG_ERR)
++ {
++ csa->stability_failures++;
++ if (csa->stability_failures >= parm->stability_lmt){
++ ret = GLP_EINSTAB;
++ goto done;
++ }
++ if (parm->msg_lev >= GLP_MSG_ERR)
+ xprintf("Warning: numerical instability (dual simplex, p"
+ "hase %s)\n", csa->phase == 1 ? "I" : "II");
+ #if 1
+@@ -3023,6 +3034,10 @@ loop: /* main loop starts here */
+ /* accuracy check based on the pivot element */
+ { double piv1 = csa->tcol_vec[csa->p]; /* more accurate */
+ double piv2 = csa->trow_vec[csa->q]; /* less accurate */
++ if(piv1 == 0.0){
++ ret = GLP_EFAIL;
++ goto done;
++ }
+ xassert(piv1 != 0.0);
+ if (fabs(piv1 - piv2) > 1e-8 * (1.0 + fabs(piv1)) ||
+ !(piv1 > 0.0 && piv2 > 0.0 || piv1 < 0.0 && piv2 < 0.0))
+--
+2.26.2
projects / cvc5.git / commitdiff