glpmpl03.c

著名的大规模线性规划求解器源码GLPK.C语言版本,可以修剪.内有详细帮助文档.

ELEMSET *set_diff
(     MPL *mpl,
      ELEMSET *X,             /* destroyed */
      ELEMSET *Y              /* destroyed */
)
{     ELEMSET *Z;
      MEMBER *memb;
      xassert(X != NULL);
      xassert(X->type == A_NONE);
      xassert(X->dim > 0);
      xassert(Y != NULL);
      xassert(Y->type == A_NONE);
      xassert(Y->dim > 0);
      xassert(X->dim == Y->dim);
      Z = create_elemset(mpl, X->dim);
      for (memb = X->head; memb != NULL; memb = memb->next)
      {  if (find_tuple(mpl, Y, memb->tuple) == NULL)
            add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple));
      }
      delete_elemset(mpl, X);
      delete_elemset(mpl, Y);
      return Z;
}

/*----------------------------------------------------------------------
-- set_symdiff - symmetric difference between two elemental sets.
--
-- This routine computes the symmetric difference:
--
--    X (+) Y = (X \ Y) U (Y \ X),
--
-- where X and Y are given elemental sets (destroyed on exit). */

ELEMSET *set_symdiff
(     MPL *mpl,
      ELEMSET *X,             /* destroyed */
      ELEMSET *Y              /* destroyed */
)
{     ELEMSET *Z;
      MEMBER *memb;
      xassert(X != NULL);
      xassert(X->type == A_NONE);
      xassert(X->dim > 0);
      xassert(Y != NULL);
      xassert(Y->type == A_NONE);
      xassert(Y->dim > 0);
      xassert(X->dim == Y->dim);
      /* Z := X \ Y */
      Z = create_elemset(mpl, X->dim);
      for (memb = X->head; memb != NULL; memb = memb->next)
      {  if (find_tuple(mpl, Y, memb->tuple) == NULL)
            add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple));
      }
      /* Z := Z U (Y \ X) */
      for (memb = Y->head; memb != NULL; memb = memb->next)
      {  if (find_tuple(mpl, X, memb->tuple) == NULL)
            add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple));
      }
      delete_elemset(mpl, X);
      delete_elemset(mpl, Y);
      return Z;
}

/*----------------------------------------------------------------------
-- set_inter - intersection of two elemental sets.
--
-- This routine computes the intersection:
--
--    X ^ Y = { j | (j in X) and (j in Y) },
--
-- where X and Y are given elemental sets (destroyed on exit). */

ELEMSET *set_inter
(     MPL *mpl,
      ELEMSET *X,             /* destroyed */
      ELEMSET *Y              /* destroyed */
)
{     ELEMSET *Z;
      MEMBER *memb;
      xassert(X != NULL);
      xassert(X->type == A_NONE);
      xassert(X->dim > 0);
      xassert(Y != NULL);
      xassert(Y->type == A_NONE);
      xassert(Y->dim > 0);
      xassert(X->dim == Y->dim);
      Z = create_elemset(mpl, X->dim);
      for (memb = X->head; memb != NULL; memb = memb->next)
      {  if (find_tuple(mpl, Y, memb->tuple) != NULL)
            add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple));
      }
      delete_elemset(mpl, X);
      delete_elemset(mpl, Y);
      return Z;
}

/*----------------------------------------------------------------------
-- set_cross - cross (Cartesian) product of two elemental sets.
--
-- This routine computes the cross (Cartesian) product:
--
--    X x Y = { (i,j) | (i in X) and (j in Y) },
--
-- where X and Y are given elemental sets (destroyed on exit). */

ELEMSET *set_cross
(     MPL *mpl,
      ELEMSET *X,             /* destroyed */
      ELEMSET *Y              /* destroyed */
)
{     ELEMSET *Z;
      MEMBER *memx, *memy;
      TUPLE *tuple, *temp;
      xassert(X != NULL);
      xassert(X->type == A_NONE);
      xassert(X->dim > 0);
      xassert(Y != NULL);
      xassert(Y->type == A_NONE);
      xassert(Y->dim > 0);
      Z = create_elemset(mpl, X->dim + Y->dim);
      for (memx = X->head; memx != NULL; memx = memx->next)
      {  for (memy = Y->head; memy != NULL; memy = memy->next)
         {  tuple = copy_tuple(mpl, memx->tuple);
            for (temp = memy->tuple; temp != NULL; temp = temp->next)
               tuple = expand_tuple(mpl, tuple, copy_symbol(mpl,
                  temp->sym));
            add_tuple(mpl, Z, tuple);
         }
      }
      delete_elemset(mpl, X);
      delete_elemset(mpl, Y);
      return Z;
}

/**********************************************************************/
/* * *                    ELEMENTAL VARIABLES                     * * */
/**********************************************************************/

/* (there are no specific routines for elemental variables) */

/**********************************************************************/
/* * *                        LINEAR FORMS                        * * */
/**********************************************************************/

/*----------------------------------------------------------------------
-- constant_term - create constant term.
--
-- This routine creates the linear form, which is a constant term. */

FORMULA *constant_term(MPL *mpl, double coef)
{     FORMULA *form;
      if (coef == 0.0)
         form = NULL;
      else
      {  form = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
         form->coef = coef;
         form->var = NULL;
         form->next = NULL;
      }
      return form;
}

/*----------------------------------------------------------------------
-- single_variable - create single variable.
--
-- This routine creates the linear form, which is a single elemental
-- variable. */

FORMULA *single_variable
(     MPL *mpl,
      ELEMVAR *var            /* referenced */
)
{     FORMULA *form;
      xassert(var != NULL);
      form = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
      form->coef = 1.0;
      form->var = var;
      form->next = NULL;
      return form;
}

/*----------------------------------------------------------------------
-- copy_formula - make copy of linear form.
--
-- This routine returns an exact copy of linear form. */

FORMULA *copy_formula
(     MPL *mpl,
      FORMULA *form           /* not changed */
)
{     FORMULA *head, *tail;
      if (form == NULL)
         head = NULL;
      else
      {  head = tail = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
         for (; form != NULL; form = form->next)
         {  tail->coef = form->coef;
            tail->var = form->var;
            if (form->next != NULL)
tail = (tail->next = dmp_get_atom(mpl->formulae, sizeof(FORMULA)));
         }
         tail->next = NULL;
      }
      return head;
}

/*----------------------------------------------------------------------
-- delete_formula - delete linear form.
--
-- This routine deletes specified linear form. */

void delete_formula
(     MPL *mpl,
      FORMULA *form           /* destroyed */
)
{     FORMULA *temp;
      while (form != NULL)
      {  temp = form;
         form = form->next;
         dmp_free_atom(mpl->formulae, temp, sizeof(FORMULA));
      }
      return;
}

/*----------------------------------------------------------------------
-- linear_comb - linear combination of two linear forms.
--
-- This routine computes the linear combination:
--
--    a * fx + b * fy,
--
-- where a and b are numeric coefficients, fx and fy are linear forms
-- (destroyed on exit). */

FORMULA *linear_comb
(     MPL *mpl,
      double a, FORMULA *fx,  /* destroyed */
      double b, FORMULA *fy   /* destroyed */
)
{     FORMULA *form = NULL, *term, *temp;
      double c0 = 0.0;
      for (term = fx; term != NULL; term = term->next)
      {  if (term->var == NULL)
            c0 = fp_add(mpl, c0, fp_mul(mpl, a, term->coef));
         else
            term->var->temp =
               fp_add(mpl, term->var->temp, fp_mul(mpl, a, term->coef));
      }
      for (term = fy; term != NULL; term = term->next)
      {  if (term->var == NULL)
            c0 = fp_add(mpl, c0, fp_mul(mpl, b, term->coef));
         else
            term->var->temp =
               fp_add(mpl, term->var->temp, fp_mul(mpl, b, term->coef));
      }
      for (term = fx; term != NULL; term = term->next)
      {  if (term->var != NULL && term->var->temp != 0.0)
         {  temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
            temp->coef = term->var->temp, temp->var = term->var;
            temp->next = form, form = temp;
            term->var->temp = 0.0;
         }
      }
      for (term = fy; term != NULL; term = term->next)
      {  if (term->var != NULL && term->var->temp != 0.0)
         {  temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
            temp->coef = term->var->temp, temp->var = term->var;
            temp->next = form, form = temp;
            term->var->temp = 0.0;
         }
      }
      if (c0 != 0.0)
      {  temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA));
         temp->coef = c0, temp->var = NULL;
         temp->next = form, form = temp;
      }
      delete_formula(mpl, fx);
      delete_formula(mpl, fy);
      return form;
}

/*----------------------------------------------------------------------
-- remove_constant - remove constant term from linear form.
--
-- This routine removes constant term from linear form and stores its
-- value to given location. */

FORMULA *remove_constant
(     MPL *mpl,
      FORMULA *form,          /* destroyed */
      double *coef            /* modified */
)
{     FORMULA *head = NULL, *temp;
      *coef = 0.0;
      while (form != NULL)
      {  temp = form;
         form = form->next;
         if (temp->var == NULL)
         {  /* constant term */
            *coef = fp_add(mpl, *coef, temp->coef);
            dmp_free_atom(mpl->formulae, temp, sizeof(FORMULA));
         }
         else
         {  /* linear term */
            temp->next = head;
            head = temp;
         }
      }
      return head;
}

/*----------------------------------------------------------------------
-- reduce_terms - reduce identical terms in linear form.
--
-- This routine reduces identical terms in specified linear form. */

FORMULA *reduce_terms
(     MPL *mpl,
      FORMULA *form           /* destroyed */
)
{     FORMULA *term, *next_term;
      double c0 = 0.0;
      for (term = form; term != NULL; term = term->next)
      {  if (term->var == NULL)
            c0 = fp_add(mpl, c0, term->coef);
         else
            term->var->temp = fp_add(mpl, term->var->temp, term->coef);
      }
      next_term = form, form = NULL;
      for (term = next_term; term != NULL; term = next_term)
      {  next_term = term->next;
         if (term->var == NULL && c0 != 0.0)
         {  term->coef = c0, c0 = 0.0;
            term->next = form, form = term;
         }
         else if (term->var != NULL && term->var->temp != 0.0)
         {  term->coef = term->var->temp, term->var->temp = 0.0;
            term->next = form, form = term;
         }
         else
            dmp_free_atom(mpl->formulae, term, sizeof(FORMULA));
      }
      return form;
}

/**********************************************************************/
/* * *                   ELEMENTAL CONSTRAINTS                    * * */
/**********************************************************************/

/* (there are no specific routines for elemental constraints) */

/**********************************************************************/
/* * *                       GENERIC VALUES                       * * */
/**********************************************************************/

/*----------------------------------------------------------------------
-- delete_value - delete generic value.
--
-- This routine deletes specified generic value.
--
-- NOTE: The generic value to be deleted must be valid. */

void delete_value
(     MPL *mpl,
      int type,
      VALUE *value            /* content destroyed */
)
{     xassert(value != NULL);
      switch (type)
      {  case A_NONE:
            value->none = NULL;
            break;
         case A_NUMERIC:
            value->num = 0.0;
            break;
         case A_SYMBOLIC:
            delete_symbol(mpl, value->sym), value->sym = NULL;
            break;
         case A_LOGICAL:
            value->bit = 0;
            break;
         case A_TUPLE:
            delete_tuple(mpl, value->tuple), value->tuple = NULL;
            break;
         case A_ELEMSET:
            delete_elemset(mpl, value->set), value->set = NULL;
            break;
         case A_ELEMVAR:
            value->var = NULL;
            break;
         case A_FORMULA:
            delete_formula(mpl, value->form), value->form = NULL;
            break;
         case A_ELEMCON:
            value->con = NULL;