--- a/NEWS
+++ b/NEWS
@@ -187,6 +187,11 @@
       helpdlg    listdlg   questdlg
       inputdlg   msgbox    warndlg
 
+ ** The glpk function has been modified to reflect changes in the GLPK
+    library. The "round" and "itcnt" options have been removed. The
+    "relax" option has been replaced by the "rtest" option. The numeric
+    values of error codes and of some options have also changed.
+
  ** Other new functions added in 3.8.0:
                                                   
       atan2d                     erfi             jit_startcnt

--- a/configure.ac
+++ b/configure.ac
@@ -869,7 +869,7 @@
 LIBS="$Z_LDFLAGS $Z_LIBS $LIBS"
 OCTAVE_CHECK_LIB(glpk, GLPK,
   [GLPK library not found.  The glpk function for solving linear programs will be disabled.],
-  [glpk/glpk.h glpk.h], [_glp_lpx_simplex])
+  [glpk/glpk.h glpk.h], [glp_simplex])
 LIBS="$save_LIBS"
 CPPFLAGS="$save_CPPFLAGS"
 

--- a/libinterp/dldfcn/__glpk__.cc
+++ b/libinterp/dldfcn/__glpk__.cc
@@ -1,6 +1,7 @@
 /*
 
 Copyright (C) 2005-2012 Nicolo' Giorgetti
+Copyright (C) 2013 Sébastien Villemot <sebastien@debian.org>
 
 This file is part of Octave.
 
@@ -46,191 +47,86 @@
 #else
 #include <glpk.h>
 #endif
-
-#if 0
-#ifdef GLPK_PRE_4_14
-
-#ifndef _GLPLIB_H
-#include <glplib.h>
-#endif
-#ifndef lib_set_fault_hook
-#define lib_set_fault_hook lib_fault_hook
-#endif
-#ifndef lib_set_print_hook
-#define lib_set_print_hook lib_print_hook
-#endif
-
-#else
-
-void _glp_lib_print_hook (int (*func)(void *info, char *buf), void *info);
-void _glp_lib_fault_hook (int (*func)(void *info, char *buf), void *info);
-
-#endif
-#endif
 }
 
-#define NIntP 17
-#define NRealP 10
-
-int lpxIntParam[NIntP] = {
-  0,
-  1,
-  0,
-  1,
-  0,
-  -1,
-  0,
-  200,
-  1,
-  2,
-  0,
-  1,
-  0,
-  0,
-  2,
-  2,
-  1
-};
-
-int IParam[NIntP] = {
-  LPX_K_MSGLEV,
-  LPX_K_SCALE,
-  LPX_K_DUAL,
-  LPX_K_PRICE,
-  LPX_K_ROUND,
-  LPX_K_ITLIM,
-  LPX_K_ITCNT,
-  LPX_K_OUTFRQ,
-  LPX_K_MPSINFO,
-  LPX_K_MPSOBJ,
-  LPX_K_MPSORIG,
-  LPX_K_MPSWIDE,
-  LPX_K_MPSFREE,
-  LPX_K_MPSSKIP,
-  LPX_K_BRANCH,
-  LPX_K_BTRACK,
-  LPX_K_PRESOL
-};
-
-
-double lpxRealParam[NRealP] = {
-  0.07,
-  1e-7,
-  1e-7,
-  1e-9,
-  -std::numeric_limits<double>::max (),
-  std::numeric_limits<double>::max (),
-  -1.0,
-  0.0,
-  1e-6,
-  1e-7
-};
-
-int RParam[NRealP] = {
-  LPX_K_RELAX,
-  LPX_K_TOLBND,
-  LPX_K_TOLDJ,
-  LPX_K_TOLPIV,
-  LPX_K_OBJLL,
-  LPX_K_OBJUL,
-  LPX_K_TMLIM,
-  LPX_K_OUTDLY,
-  LPX_K_TOLINT,
-  LPX_K_TOLOBJ
-};
+typedef struct
+{
+  int msglev;
+  int dual;
+  int price;
+  int itlim;
+  int outfrq;
+  int branch;
+  int btrack;
+  int presol;
+  int rtest;
+  int tmlim;
+  int outdly;
+  double tolbnd;
+  double toldj;
+  double tolpiv;
+  double objll;
+  double objul;
+  double tolint;
+  double tolobj;
+} control_params;
 
 static jmp_buf mark;  //-- Address for long jump to jump to
 
-#if 0
-int
-glpk_fault_hook (void * /* info */, char *msg)
-{
-  error ("CRITICAL ERROR in GLPK: %s", msg);
-  longjmp (mark, -1);
-}
-
-int
-glpk_print_hook (void * /* info */, char *msg)
-{
-  message (0, "%s", msg);
-  return 1;
-}
-#endif
-
 int
 glpk (int sense, int n, int m, double *c, int nz, int *rn, int *cn,
       double *a, double *b, char *ctype, int *freeLB, double *lb,
       int *freeUB, double *ub, int *vartype, int isMIP, int lpsolver,
-      int save_pb, double *xmin, double *fmin, double *status,
-      double *lambda, double *redcosts, double *time, double *mem)
+      int save_pb, int scale, const control_params *par,
+      double *xmin, double *fmin, int *status,
+      double *lambda, double *redcosts, double *time)
 {
-  int errnum;
   int typx = 0;
-  int method;
+  int errnum = 0;
 
   clock_t t_start = clock ();
 
-#if 0
-#ifdef GLPK_PRE_4_14
-  lib_set_fault_hook (0, glpk_fault_hook);
-#else
-  _glp_lib_fault_hook (glpk_fault_hook, 0);
-#endif
-
-  if (lpxIntParam[0] > 1)
-#ifdef GLPK_PRE_4_14
-    lib_set_print_hook (0, glpk_print_hook);
-#else
-    _glp_lib_print_hook (glpk_print_hook, 0);
-#endif
-#endif
-
-  LPX *lp = lpx_create_prob ();
-
+  glp_prob *lp = glp_create_prob ();
 
   //-- Set the sense of optimization
   if (sense == 1)
-    lpx_set_obj_dir (lp, LPX_MIN);
+    glp_set_obj_dir (lp, GLP_MIN);
   else
-    lpx_set_obj_dir (lp, LPX_MAX);
+    glp_set_obj_dir (lp, GLP_MAX);
 
-  //-- If the problem has integer structural variables switch to MIP
-  if (isMIP)
-    lpx_set_class (lp, LPX_MIP);
-
-  lpx_add_cols (lp, n);
+  glp_add_cols (lp, n);
   for (int i = 0; i < n; i++)
     {
       //-- Define type of the structural variables
       if (! freeLB[i] && ! freeUB[i])
         {
           if (lb[i] != ub[i])
-            lpx_set_col_bnds (lp, i+1, LPX_DB, lb[i], ub[i]);
+            glp_set_col_bnds (lp, i+1, GLP_DB, lb[i], ub[i]);
           else
-            lpx_set_col_bnds (lp, i+1, LPX_FX, lb[i], ub[i]);
+            glp_set_col_bnds (lp, i+1, GLP_FX, lb[i], ub[i]);
         }
       else
         {
           if (! freeLB[i] && freeUB[i])
-            lpx_set_col_bnds (lp, i+1, LPX_LO, lb[i], ub[i]);
+            glp_set_col_bnds (lp, i+1, GLP_LO, lb[i], ub[i]);
           else
             {
               if (freeLB[i] && ! freeUB[i])
-                lpx_set_col_bnds (lp, i+1, LPX_UP, lb[i], ub[i]);
+                glp_set_col_bnds (lp, i+1, GLP_UP, lb[i], ub[i]);
               else
-                lpx_set_col_bnds (lp, i+1, LPX_FR, lb[i], ub[i]);
+                glp_set_col_bnds (lp, i+1, GLP_FR, lb[i], ub[i]);
             }
         }
 
       // -- Set the objective coefficient of the corresponding
       // -- structural variable. No constant term is assumed.
-      lpx_set_obj_coef(lp,i+1,c[i]);
+      glp_set_obj_coef(lp,i+1,c[i]);
 
       if (isMIP)
-        lpx_set_col_kind (lp, i+1, vartype[i]);
+        glp_set_col_kind (lp, i+1, vartype[i]);
     }
 
-  lpx_add_rows (lp, m);
+  glp_add_rows (lp, m);
 
   for (int i = 0; i < m; i++)
     {
@@ -245,124 +141,123 @@
       switch (ctype[i])
         {
         case 'F':
-          typx = LPX_FR;
+          typx = GLP_FR;
           break;
 
         case 'U':
-          typx = LPX_UP;
+          typx = GLP_UP;
           break;
 
         case 'L':
-          typx = LPX_LO;
+          typx = GLP_LO;
           break;
 
         case 'S':
-          typx = LPX_FX;
+          typx = GLP_FX;
           break;
 
         case 'D':
-          typx = LPX_DB;
+          typx = GLP_DB;
           break;
         }
 
-      lpx_set_row_bnds (lp, i+1, typx, b[i], b[i]);
+      glp_set_row_bnds (lp, i+1, typx, b[i], b[i]);
 
     }
 
-  lpx_load_matrix (lp, nz, rn, cn, a);
+  glp_load_matrix (lp, nz, rn, cn, a);
 
   if (save_pb)
     {
       static char tmp[] = "outpb.lp";
-      if (lpx_write_cpxlp (lp, tmp) != 0)
+      if (glp_write_lp (lp, NULL, tmp) != 0)
         {
           error ("__glpk__: unable to write problem");
           longjmp (mark, -1);
         }
     }
 
-  //-- scale the problem data (if required)
-  //-- if (scale && (!presol || method == 1)) lpx_scale_prob (lp);
-  //-- LPX_K_SCALE=IParam[1]  LPX_K_PRESOL=IParam[16]
-  if (lpxIntParam[1] && (! lpxIntParam[16] || lpsolver != 1))
-    lpx_scale_prob (lp);
+  //-- scale the problem data
+  if (!par->presol || lpsolver != 1)
+    glp_scale_prob (lp, scale);
 
   //-- build advanced initial basis (if required)
-  if (lpsolver == 1 && ! lpxIntParam[16])
-    lpx_adv_basis (lp);
-
-  for (int i = 0; i < NIntP; i++)
-    lpx_set_int_parm (lp, IParam[i], lpxIntParam[i]);
+  if (lpsolver == 1 && !par->presol)
+    glp_adv_basis (lp, 0);
 
-  for (int i = 0; i < NRealP; i++)
-    lpx_set_real_parm (lp, RParam[i], lpxRealParam[i]);
-
-  if (lpsolver == 1)
-    method = 'S';
-  else
-    method = 'T';
-
-  switch (method)
+  /* For MIP problems without a presolver, a first pass with glp_simplex
+     is required */
+  if ((!isMIP && lpsolver == 1)
+      || (isMIP && !par->presol))
     {
-    case 'S':
-      {
-        if (isMIP)
-          {
-            method = 'I';
-            errnum = lpx_simplex (lp);
-            errnum = lpx_integer (lp);
-          }
-        else
-          errnum = lpx_simplex (lp);
-      }
-     break;
-
-    case 'T':
-      errnum = lpx_interior (lp);
-      break;
-
-    default:
-      break;
-#if 0
-#ifdef GLPK_PRE_4_14
-      insist (method != method);
-#else
-      static char tmp[] = "method != method";
-      glpk_fault_hook (0, tmp);
-#endif
-#endif
+      glp_smcp smcp;
+      glp_init_smcp (&smcp);
+      smcp.msg_lev = par->msglev;
+      smcp.meth = par->dual;
+      smcp.pricing = par->price;
+      smcp.r_test = par->rtest;
+      smcp.tol_bnd = par->tolbnd;
+      smcp.tol_dj = par->toldj;
+      smcp.tol_piv = par->tolpiv;
+      smcp.obj_ll = par->objll;
+      smcp.obj_ul = par->objul;
+      smcp.it_lim = par->itlim;
+      smcp.tm_lim = par->tmlim;
+      smcp.out_frq = par->outfrq;
+      smcp.out_dly = par->outdly;
+      smcp.presolve = par->presol;
+      errnum = glp_simplex (lp, &smcp);
     }
 
-  /*  errnum assumes the following results:
-      errnum = 0 <=> No errors
-      errnum = 1 <=> Iteration limit exceeded.
-      errnum = 2 <=> Numerical problems with basis matrix.
-  */
-  if (errnum == LPX_E_OK)
+  if (isMIP)
+    {
+      glp_iocp iocp;
+      glp_init_iocp (&iocp);
+      iocp.msg_lev = par->msglev;
+      iocp.br_tech = par->branch;
+      iocp.bt_tech = par->btrack;
+      iocp.tol_int = par->tolint;
+      iocp.tol_obj = par->tolobj;
+      iocp.tm_lim = par->tmlim;
+      iocp.out_frq = par->outfrq;
+      iocp.out_dly = par->outdly;
+      iocp.presolve = par->presol;
+      errnum = glp_intopt (lp, &iocp);
+    }
+
+  if (!isMIP && lpsolver == 2)
+    {
+      glp_iptcp iptcp;
+      glp_init_iptcp (&iptcp);
+      iptcp.msg_lev = par->msglev;
+      errnum = glp_interior (lp, &iptcp);
+    }
+
+  if (errnum == 0)
     {
       if (isMIP)
         {
-          *status = lpx_mip_status (lp);
-          *fmin = lpx_mip_obj_val (lp);
+          *status = glp_mip_status (lp);
+          *fmin = glp_mip_obj_val (lp);
         }
       else
         {
           if (lpsolver == 1)
             {
-              *status = lpx_get_status (lp);
-              *fmin = lpx_get_obj_val (lp);
+              *status = glp_get_status (lp);
+              *fmin = glp_get_obj_val (lp);
             }
           else
             {
-              *status = lpx_ipt_status (lp);
-              *fmin = lpx_ipt_obj_val (lp);
+              *status = glp_ipt_status (lp);
+              *fmin = glp_ipt_obj_val (lp);
             }
         }
 
       if (isMIP)
         {
           for (int i = 0; i < n; i++)
-            xmin[i] = lpx_mip_col_val (lp, i+1);
+            xmin[i] = glp_mip_col_val (lp, i+1);
         }
       else
         {
@@ -370,52 +265,41 @@
           for (int i = 0; i < n; i++)
             {
               if (lpsolver == 1)
-                xmin[i] = lpx_get_col_prim (lp, i+1);
+                xmin[i] = glp_get_col_prim (lp, i+1);
               else
-                xmin[i] = lpx_ipt_col_prim (lp, i+1);
+                xmin[i] = glp_ipt_col_prim (lp, i+1);
             }
 
           /* Dual values */
           for (int i = 0; i < m; i++)
             {
               if (lpsolver == 1)
-                lambda[i] = lpx_get_row_dual (lp, i+1);
+                lambda[i] = glp_get_row_dual (lp, i+1);
               else
-                lambda[i] = lpx_ipt_row_dual (lp, i+1);
+                lambda[i] = glp_ipt_row_dual (lp, i+1);
             }
 
           /* Reduced costs */
-          for (int i = 0; i < lpx_get_num_cols (lp); i++)
+          for (int i = 0; i < glp_get_num_cols (lp); i++)
             {
               if (lpsolver == 1)
-                redcosts[i] = lpx_get_col_dual (lp, i+1);
+                redcosts[i] = glp_get_col_dual (lp, i+1);
               else
-                redcosts[i] = lpx_ipt_col_dual (lp, i+1);
+                redcosts[i] = glp_ipt_col_dual (lp, i+1);
             }
         }
 
       *time = (clock () - t_start) / CLOCKS_PER_SEC;
-
-#ifdef GLPK_PRE_4_14
-      *mem = (lib_env_ptr () -> mem_tpeak);
-#else
-      *mem = 0;
-#endif
-
-      lpx_delete_prob (lp);
-      return 0;
     }
 
-   lpx_delete_prob (lp);
-
-   *status = errnum;
+   glp_delete_prob (lp);
 
    return errnum;
 }
 
 #endif
 
-#define OCTAVE_GLPK_GET_REAL_PARAM(NAME, IDX) \
+#define OCTAVE_GLPK_GET_REAL_PARAM(NAME, VAL) \
   do \
     { \
       octave_value tmp = PARAM.getfield (NAME); \
@@ -424,7 +308,7 @@
         { \
           if (! tmp.is_empty ()) \
             { \
-              lpxRealParam[IDX] = tmp.scalar_value (); \
+              VAL = tmp.scalar_value (); \
  \
               if (error_state) \
                 { \
@@ -659,10 +543,10 @@
       if (VTYPE(i,0) == 'I')
         {
           isMIP = 1;
-          vartype(i) = LPX_IV;
+          vartype(i) = GLP_IV;
         }
       else
-        vartype(i) = LPX_CV;
+        vartype(i) = GLP_CV;
     }
 
   //-- 8th Input. Sense of optimization.
@@ -689,78 +573,78 @@
       return retval;
     }
 
+  control_params par;
+
   //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   //-- Integer parameters
   //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
   //-- Level of messages output by the solver
-  OCTAVE_GLPK_GET_INT_PARAM ("msglev", lpxIntParam[0]);
-  if (lpxIntParam[0] < 0 || lpxIntParam[0] > 3)
+  par.msglev = 1;
+  OCTAVE_GLPK_GET_INT_PARAM ("msglev", par.msglev);
+  if (par.msglev < 0 || par.msglev > 3)
     {
-      error ("__glpk__: PARAM.msglev must be 0 (no output [default]) or 1 (error messages only) or 2 (normal output) or 3 (full output)");
+      error ("__glpk__: PARAM.msglev must be 0 (no output) or 1 (error and warning messages only [default]) or 2 (normal output) or 3 (full output)");
       return retval;
     }
 
   //-- scaling option
-  OCTAVE_GLPK_GET_INT_PARAM ("scale", lpxIntParam[1]);
-  if (lpxIntParam[1] < 0 || lpxIntParam[1] > 2)
+  volatile int scale = 16;
+  OCTAVE_GLPK_GET_INT_PARAM ("scale", scale);
+  if (scale < 0 || scale > 128)
     {
-      error ("__glpk__: PARAM.scale must be 0 (no scaling) or 1 (equilibration scaling [default]) or 2 (geometric mean scaling)");
+      error ("__glpk__: PARAM.scale must either be 128 (automatic selection of scaling options), or a bitwise or of: 1 (geometric mean scaling), 16 (equilibration scaling), 32 (round scale factors to power of two), 64 (skip if problem is well scaled");
       return retval;
     }
 
-  //-- Dual dimplex option
-  OCTAVE_GLPK_GET_INT_PARAM ("dual", lpxIntParam[2]);
-  if (lpxIntParam[2] < 0 || lpxIntParam[2] > 1)
+  //-- Dual simplex option
+  par.dual = 1;
+  OCTAVE_GLPK_GET_INT_PARAM ("dual", par.dual);
+  if (par.dual < 1 || par.dual > 3)
     {
-      error ("__glpk__: PARAM.dual must be 0 (do NOT use dual simplex [default]) or 1 (use dual simplex)");
+      error ("__glpk__: PARAM.dual must be 1 (use two-phase primal simplex [default]) or 2 (use two-phase dual simplex) or 3 (use two-phase dual simplex, and if it fails, switch to the primal simplex)");
       return retval;
     }
 
   //-- Pricing option
-  OCTAVE_GLPK_GET_INT_PARAM ("price", lpxIntParam[3]);
-  if (lpxIntParam[3] < 0 || lpxIntParam[3] > 1)
+  par.price = 34;
+  OCTAVE_GLPK_GET_INT_PARAM ("price", par.price);
+  if (par.price != 17 && par.price != 34)
     {
-      error ("__glpk__: PARAM.price must be 0 (textbook pricing) or 1 (steepest edge pricing [default])");
-      return retval;
-    }
-
-  //-- Solution rounding option
-  OCTAVE_GLPK_GET_INT_PARAM ("round", lpxIntParam[4]);
-  if (lpxIntParam[4] < 0 || lpxIntParam[4] > 1)
-    {
-      error ("__glpk__: PARAM.round must be 0 (report all primal and dual values [default]) or 1 (replace tiny primal and dual values by exact zero)");
+      error ("__glpk__: PARAM.price must be 17 (textbook pricing) or 34 (steepest edge pricing [default])");
       return retval;
     }
 
   //-- Simplex iterations limit
-  OCTAVE_GLPK_GET_INT_PARAM ("itlim", lpxIntParam[5]);
-
-  //-- Simplex iterations count
-  OCTAVE_GLPK_GET_INT_PARAM ("itcnt", lpxIntParam[6]);
+  par.itlim = INT_MAX;
+  OCTAVE_GLPK_GET_INT_PARAM ("itlim", par.itlim);
 
   //-- Output frequency, in iterations
-  OCTAVE_GLPK_GET_INT_PARAM ("outfrq", lpxIntParam[7]);
+  par.outfrq = 200;
+  OCTAVE_GLPK_GET_INT_PARAM ("outfrq", par.outfrq);
 
   //-- Branching heuristic option
-  OCTAVE_GLPK_GET_INT_PARAM ("branch", lpxIntParam[14]);
-  if (lpxIntParam[14] < 0 || lpxIntParam[14] > 2)
+  par.branch = 4;
+  OCTAVE_GLPK_GET_INT_PARAM ("branch", par.branch);
+  if (par.branch < 1 || par.branch > 5)
     {
-      error ("__glpk__: PARAM.branch must be (MIP only) 0 (branch on first variable) or 1 (branch on last variable) or 2 (branch using a heuristic by Driebeck and Tomlin [default]");
+      error ("__glpk__: PARAM.branch must be 1 (first fractional variable) or 2 (last fractional variable) or 3 (most fractional variable) or 4 (heuristic by Driebeck and Tomlin [default]) or 5 (hybrid pseudocost heuristic)");
       return retval;
     }
 
   //-- Backtracking heuristic option
-  OCTAVE_GLPK_GET_INT_PARAM ("btrack", lpxIntParam[15]);
-  if (lpxIntParam[15] < 0 || lpxIntParam[15] > 2)
+  par.btrack = 4;
+  OCTAVE_GLPK_GET_INT_PARAM ("btrack", par.btrack);
+  if (par.btrack < 1 || par.btrack > 4)
     {
-      error ("__glpk__: PARAM.btrack must be (MIP only) 0 (depth first search) or 1 (breadth first search) or 2 (backtrack using the best projection heuristic [default]");
+      error ("__glpk__: PARAM.btrack must be 1 (depth first search) or 2 (breadth first search) or 3 (best local bound) or 4 (best projection heuristic [default]");
       return retval;
     }
 
   //-- Presolver option
-  OCTAVE_GLPK_GET_INT_PARAM ("presol", lpxIntParam[16]);
-  if (lpxIntParam[16] < 0 || lpxIntParam[16] > 1)
+  par.presol = 1;
+  OCTAVE_GLPK_GET_INT_PARAM ("presol", par.presol);
+  if (par.presol < 0 || par.presol > 1)
     {
       error ("__glpk__: PARAM.presol must be 0 (do NOT use LP presolver) or 1 (use LP presolver [default])");
       return retval;
@@ -775,6 +659,21 @@
       return retval;
     }
 
+  //-- Ratio test option
+  par.rtest = 34;
+  OCTAVE_GLPK_GET_INT_PARAM ("rtest", par.rtest);
+  if (par.rtest != 17 && par.rtest != 34)
+    {
+      error ("__glpk__: PARAM.rtest must be 17 (standard ratio test) or 34 (Harris' two-pass ratio test [default])");
+      return retval;
+    }
+
+  par.tmlim = INT_MAX;
+  OCTAVE_GLPK_GET_INT_PARAM ("tmlim", par.tmlim);
+
+  par.outdly = 0;
+  OCTAVE_GLPK_GET_INT_PARAM ("outdly", par.outdly);
+
   //-- Save option
   volatile int save_pb = 0;
   OCTAVE_GLPK_GET_INT_PARAM ("save", save_pb);
@@ -784,51 +683,50 @@
   //-- Real parameters
   //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
-  //-- Ratio test option
-  OCTAVE_GLPK_GET_REAL_PARAM ("relax", 0);
-
   //-- Relative tolerance used to check if the current basic solution
   //-- is primal feasible
-  OCTAVE_GLPK_GET_REAL_PARAM ("tolbnd", 1);
+  par.tolbnd = 1e-7;
+  OCTAVE_GLPK_GET_REAL_PARAM ("tolbnd", par.tolbnd);
 
   //-- Absolute tolerance used to check if the current basic solution
   //-- is dual feasible
-  OCTAVE_GLPK_GET_REAL_PARAM ("toldj", 2);
+  par.toldj = 1e-7;
+  OCTAVE_GLPK_GET_REAL_PARAM ("toldj", par.toldj);
 
   //-- Relative tolerance used to choose eligible pivotal elements of
   //--  the simplex table in the ratio test
-  OCTAVE_GLPK_GET_REAL_PARAM ("tolpiv", 3);
+  par.tolpiv = 1e-10;
+  OCTAVE_GLPK_GET_REAL_PARAM ("tolpiv", par.tolpiv);
 
-  OCTAVE_GLPK_GET_REAL_PARAM ("objll", 4);
-
-  OCTAVE_GLPK_GET_REAL_PARAM ("objul", 5);
+  par.objll = -std::numeric_limits<double>::max ();
+  OCTAVE_GLPK_GET_REAL_PARAM ("objll", par.objll);
 
-  OCTAVE_GLPK_GET_REAL_PARAM ("tmlim", 6);
-
-  OCTAVE_GLPK_GET_REAL_PARAM ("outdly", 7);
+  par.objul = std::numeric_limits<double>::max ();
+  OCTAVE_GLPK_GET_REAL_PARAM ("objul", par.objul);
 
-  OCTAVE_GLPK_GET_REAL_PARAM ("tolint", 8);
+  par.tolint = 1e-5;
+  OCTAVE_GLPK_GET_REAL_PARAM ("tolint", par.tolint);
 
-  OCTAVE_GLPK_GET_REAL_PARAM ("tolobj", 9);
+  par.tolobj = 1e-7;
+  OCTAVE_GLPK_GET_REAL_PARAM ("tolobj", par.tolobj);
 
   //-- Assign pointers to the output parameters
   ColumnVector xmin (mrowsc, octave_NA);
   double fmin = octave_NA;
-  double status;
   ColumnVector lambda (mrowsA, octave_NA);
   ColumnVector redcosts (mrowsc, octave_NA);
   double time;
-  double mem;
+  int status, errnum = 0;
 
   int jmpret = setjmp (mark);
 
   if (jmpret == 0)
-    glpk (sense, mrowsc, mrowsA, c, nz, rn.fortran_vec (),
-          cn.fortran_vec (), a.fortran_vec (), b, ctype,
-          freeLB.fortran_vec (), lb, freeUB.fortran_vec (), ub,
-          vartype.fortran_vec (), isMIP, lpsolver, save_pb,
-          xmin.fortran_vec (), &fmin, &status, lambda.fortran_vec (),
-          redcosts.fortran_vec (), &time, &mem);
+    errnum = glpk (sense, mrowsc, mrowsA, c, nz, rn.fortran_vec (),
+                   cn.fortran_vec (), a.fortran_vec (), b, ctype,
+                   freeLB.fortran_vec (), lb, freeUB.fortran_vec (), ub,
+                   vartype.fortran_vec (), isMIP, lpsolver, save_pb, scale, &par,
+                   xmin.fortran_vec (), &fmin, &status, lambda.fortran_vec (),
+                   redcosts.fortran_vec (), &time);
 
   octave_scalar_map extra;
 
@@ -839,10 +737,10 @@
     }
 
   extra.assign ("time", time);
-  extra.assign ("mem", mem);
+  extra.assign ("status", status);
 
   retval(3) = extra;
-  retval(2) = status;
+  retval(2) = errnum;
   retval(1) = fmin;
   retval(0) = xmin;
 

--- a/scripts/optimization/glpk.m
+++ b/scripts/optimization/glpk.m
@@ -1,4 +1,5 @@
 ## Copyright (C) 2005-2012 Nicolo' Giorgetti
+## Copyright (C) 2013 Sébastien Villemot <sebastien@debian.org>
 ##
 ## This file is part of Octave.
 ##
@@ -17,7 +18,7 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn {Function File} {[@var{xopt}, @var{fmin}, @var{status}, @var{extra}] =} glpk (@var{c}, @var{A}, @var{b}, @var{lb}, @var{ub}, @var{ctype}, @var{vartype}, @var{sense}, @var{param})
+## @deftypefn {Function File} {[@var{xopt}, @var{fmin}, @var{errnum}, @var{extra}] =} glpk (@var{c}, @var{A}, @var{b}, @var{lb}, @var{ub}, @var{ctype}, @var{vartype}, @var{sense}, @var{param})
 ## Solve a linear program using the GNU @sc{glpk} library.  Given three
 ## arguments, @code{glpk} solves the following standard LP:
 ## @tex
@@ -149,112 +150,119 @@
 ## Integer parameters:
 ##
 ## @table @code
-## @item msglev (@w{@code{LPX_K_MSGLEV}}, default: 1)
+## @item msglev (default: 1)
 ## Level of messages output by solver routines:
 ##
 ## @table @asis
-## @item 0
+## @item 0 (@w{@code{GLP_MSG_OFF}})
 ## No output.
 ##
-## @item 1
-## Error messages only.
+## @item 1 (@w{@code{GLP_MSG_ERR}})
+## Error and warning messages only.
 ##
-## @item 2
+## @item 2 (@w{@code{GLP_MSG_ON}})
 ## Normal output.
 ##
-## @item 3
+## @item 3 (@w{@code{GLP_MSG_ALL}})
 ## Full output (includes informational messages).
 ## @end table
 ##
-## @item scale (@w{@code{LPX_K_SCALE}}, default: 1)
-## Scaling option:
+## @item scale (default: 16)
+## Scaling option. The values can be combined with the bitwise OR operator and
+## may be the following:
 ##
 ## @table @asis
-## @item 0
-## No scaling.
+## @item 1 (@w{@code{GLP_SF_GM}})
+## Geometric mean scaling.
 ##
-## @item 1
+## @item 16 (@w{@code{GLP_SF_EQ}})
 ## Equilibration scaling.
 ##
-## @item 2
-## Geometric mean scaling, then equilibration scaling.
+## @item 32 (@w{@code{GLP_SF_2N}})
+## Round scale factors to power of two.
+##
+## @item 64 (@w{@code{GLP_SF_SKIP}})
+## Skip if problem is well scaled.
 ## @end table
 ##
-## @item dual    (@w{@code{LPX_K_DUAL}}, default: 0)
-## Dual simplex option:
+## Alternatively, a value of 128 (@code{GLP_SF_AUTO}) may be also specified, in which case the
+## routine chooses the scaling options automatically.
+##
+## @item dual (default: 1)
+## Simplex method option:
 ##
 ## @table @asis
-## @item 0
-## Do not use the dual simplex.
+## @item 1 (@w{@code{GLP_PRIMAL}})
+## Use two-phase primal simplex.
 ##
-## @item 1
-## If initial basic solution is dual feasible, use the dual simplex.
+## @item 2 (@w{@code{GLP_DUALP}})
+## Use two-phase dual simplex, and if it fails, switch to the primal simplex.
+##
+## @item 3 (@w{@code{GLP_DUAL}})
+## Use two-phase dual simplex.
 ## @end table
 ##
-## @item price   (@w{@code{LPX_K_PRICE}}, default: 1)
+## @item price (default: 34)
 ## Pricing option (for both primal and dual simplex):
 ##
 ## @table @asis
-## @item 0
+## @item 17 (@w{@code{GLP_PT_STD}})
 ## Textbook pricing.
 ##
-## @item 1
+## @item 34 (@w{@code{GLP_PT_PSE}})
 ## Steepest edge pricing.
 ## @end table
 ##
-## @item round   (@w{@code{LPX_K_ROUND}}, default: 0)
-## Solution rounding option:
-##
-## @table @asis
-## @item 0
-## Report all primal and dual values "as is".
+## @item itlim (default: intmax)
+## Simplex iterations limit.  It is decreased by one each time when one simplex
+## iteration has been performed, and reaching zero value signals the solver to
+## stop the search.
 ##
-## @item 1
-## Replace tiny primal and dual values by exact zero.
-## @end table
-##
-## @item itlim   (@w{@code{LPX_K_ITLIM}}, default: -1)
-## Simplex iterations limit.  If this value is positive, it is decreased by
-## one each time when one simplex iteration has been performed, and
-## reaching zero value signals the solver to stop the search.  Negative
-## value means no iterations limit.
-##
-## @item itcnt (@w{@code{LPX_K_OUTFRQ}}, default: 200)
+## @item outfrq (default: 200)
 ## Output frequency, in iterations.  This parameter specifies how
 ## frequently the solver sends information about the solution to the
 ## standard output.
 ##
-## @item branch (@w{@code{LPX_K_BRANCH}}, default: 2)
-## Branching heuristic option (for MIP only):
+## @item branch (default: 4)
+## Branching technique option (for MIP only):
 ##
 ## @table @asis
-## @item 0
-## Branch on the first variable.
+## @item 1 (@w{@code{GLP_BR_FFV}})
+## First fractional variable.
+##
+## @item 2 (@w{@code{GLP_BR_LFV}})
+## Last fractional variable.
 ##
-## @item 1
-## Branch on the last variable.
+## @item 3 (@w{@code{GLP_BR_MFV}})
+## Most fractional variable.
 ##
-## @item 2
-## Branch using a heuristic by Driebeck and Tomlin.
+## @item 4 (@w{@code{GLP_BR_DTH}})
+## Heuristic by Driebeck and Tomlin.
+##
+## @item 5 (@w{@code{GLP_BR_PCH}})
+## Hybrid pseudocost heuristic.
 ## @end table
 ##
-## @item btrack (@w{@code{LPX_K_BTRACK}}, default: 2)
-## Backtracking heuristic option (for MIP only):
+## @item btrack (default: 4)
+## Backtracking technique option (for MIP only):
 ##
 ## @table @asis
-## @item 0
+## @item 1 (@w{@code{GLP_BT_DFS}})
 ## Depth first search.
 ##
-## @item 1
+## @item 2 (@w{@code{GLP_BT_BFS}})
 ## Breadth first search.
 ##
-## @item 2
-## Backtrack using the best projection heuristic.
+## @item 3 (@w{@code{GLP_BT_BLB}})
+## Best local bound.
+##
+## @item 4 (@w{@code{GLP_BT_BPH}})
+## Best projection heuristic.
 ## @end table
 ##
-## @item presol (@w{@code{LPX_K_PRESOL}}, default: 1)
-## If this flag is set, the routine lpx_simplex solves the problem using
-## the built-in LP presolver.  Otherwise the LP presolver is not used.
+## @item presol (default: 1)
+## If this flag is set, the simplex solver uses the built-in LP presolver.
+## Otherwise the LP presolver is not used.
 ##
 ## @item lpsolver (default: 1)
 ## Select which solver to use.  If the problem is a MIP problem this flag
@@ -268,6 +276,25 @@
 ## Interior point method.
 ## @end table
 ##
+## @item rtest (default: 34)
+## Ratio test technique:
+##
+## @table @asis
+## @item 17 (@w{@code{GLP_RT_STD}})
+## Standard ("textbook").
+##
+## @item 34 (@w{@code{GLP_RT_HAR}})
+## Harris' two-pass ratio test.
+## @end table
+##
+## @item tmlim (default: intmax)
+## Searching time limit, in milliseconds.
+##
+## @item outdly (default: 0)
+## Output delay, in seconds.  This parameter specifies how long the solver
+## should delay sending information about the solution to the standard
+## output.
+##
 ## @item save (default: 0)
 ## If this parameter is nonzero, save a copy of the problem in
 ## CPLEX LP format to the file @file{"outpb.lp"}.  There is currently no
@@ -277,58 +304,37 @@
 ## Real parameters:
 ##
 ## @table @code
-## @item relax (@w{@code{LPX_K_RELAX}}, default: 0.07)
-## Relaxation parameter used in the ratio test.  If it is zero, the textbook
-## ratio test is used.  If it is non-zero (should be positive), Harris'
-## two-pass ratio test is used.  In the latter case on the first pass of the
-## ratio test basic variables (in the case of primal simplex) or reduced
-## costs of non-basic variables (in the case of dual simplex) are allowed
-## to slightly violate their bounds, but not more than
-## @code{relax*tolbnd} or @code{relax*toldj (thus, @code{relax} is a
-## percentage of @code{tolbnd} or @code{toldj}}.
-##
-## @item tolbnd (@w{@code{LPX_K_TOLBND}}, default: 10e-7)
+## @item tolbnd (default: 1e-7)
 ## Relative tolerance used to check if the current basic solution is primal
 ## feasible.  It is not recommended that you change this parameter unless you
 ## have a detailed understanding of its purpose.
 ##
-## @item toldj (@w{@code{LPX_K_TOLDJ}}, default: 10e-7)
+## @item toldj (default: 1e-7)
 ## Absolute tolerance used to check if the current basic solution is dual
 ## feasible.  It is not recommended that you change this parameter unless you
 ## have a detailed understanding of its purpose.
 ##
-## @item tolpiv (@w{@code{LPX_K_TOLPIV}}, default: 10e-9)
+## @item tolpiv (default: 1e-10)
 ## Relative tolerance used to choose eligible pivotal elements of the
 ## simplex table.  It is not recommended that you change this parameter unless
 ## you have a detailed understanding of its purpose.
 ##
-## @item objll (@w{@code{LPX_K_OBJLL}}, default: -DBL_MAX)
-## Lower limit of the objective function.  If on the phase II the objective
+## @item objll (default: -DBL_MAX)
+## Lower limit of the objective function.  If the objective
 ## function reaches this limit and continues decreasing, the solver stops
 ## the search.  This parameter is used in the dual simplex method only.
 ##
-## @item objul (@w{@code{LPX_K_OBJUL}}, default: +DBL_MAX)
-## Upper limit of the objective function.  If on the phase II the objective
+## @item objul (default: +DBL_MAX)
+## Upper limit of the objective function.  If the objective
 ## function reaches this limit and continues increasing, the solver stops
 ## the search.  This parameter is used in the dual simplex only.
 ##
-## @item tmlim (@w{@code{LPX_K_TMLIM}}, default: -1.0)
-## Searching time limit, in seconds.  If this value is positive, it is
-## decreased each time when one simplex iteration has been performed by the
-## amount of time spent for the iteration, and reaching zero value signals
-## the solver to stop the search.  Negative value means no time limit.
-##
-## @item outdly (@w{@code{LPX_K_OUTDLY}}, default: 0.0)
-## Output delay, in seconds.  This parameter specifies how long the solver
-## should delay sending information about the solution to the standard
-## output.  Non-positive value means no delay.
-##
-## @item tolint (@w{@code{LPX_K_TOLINT}}, default: 10e-5)
+## @item tolint (default: 1e-5)
 ## Relative tolerance used to check if the current basic solution is integer
 ## feasible.  It is not recommended that you change this parameter unless
 ## you have a detailed understanding of its purpose.
 ##
-## @item tolobj (@w{@code{LPX_K_TOLOBJ}}, default: 10e-7)
+## @item tolobj (default: 1e-7)
 ## Relative tolerance used to check if the value of the objective function
 ## is not better than in the best known integer feasible solution.  It is
 ## not recommended that you change this parameter unless you have a
@@ -345,94 +351,69 @@
 ## @item fopt
 ## The optimum value of the objective function.
 ##
-## @item status
-## Status of the optimization.
-##
-## Simplex Method:
-##
-## @table @asis
-## @item 180 (@w{@code{LPX_OPT}})
-## Solution is optimal.
-##
-## @item 181 (@w{@code{LPX_FEAS}})
-## Solution is feasible.
-##
-## @item 182 (@w{@code{LPX_INFEAS}})
-## Solution is infeasible.
-##
-## @item 183 (@w{@code{LPX_NOFEAS}})
-## Problem has no feasible solution.
-##
-## @item 184 (@w{@code{LPX_UNBND}})
-## Problem has no unbounded solution.
-##
-## @item 185 (@w{@code{LPX_UNDEF}})
-## Solution status is undefined.
-## @end table
-##
-## Interior Point Method:
-##
-## @table @asis
-## @item 150 (@w{@code{LPX_T_UNDEF}})
-## The interior point method is undefined.
-##
-## @item 151 (@w{@code{LPX_T_OPT}})
-## The interior point method is optimal.
-## @end table
-##
-## Mixed Integer Method:
+## @item errnum
+## Error code.
 ##
 ## @table @asis
-## @item 170 (@w{@code{LPX_I_UNDEF}})
-## The status is undefined.
+## @item 0
+## No error.
 ##
-## @item 171 (@w{@code{LPX_I_OPT}})
-## The solution is integer optimal.
+## @item 1 (@w{@code{GLP_EBADB}})
+## Invalid basis.
 ##
-## @item 172 (@w{@code{LPX_I_FEAS}})
-## Solution integer feasible but its optimality has not been proven
+## @item 2 (@w{@code{GLP_ESING}})
+## Singular matrix.
 ##
-## @item 173 (@w{@code{LPX_I_NOFEAS}})
-## No integer feasible solution.
-## @end table
+## @item 3 (@w{@code{GLP_ECOND}})
+## Ill-conditioned matrix.
 ##
-## @noindent
-## If an error occurs, @var{status} will contain one of the following
-## codes:
+## @item 4 (@w{@code{GLP_EBOUND}})
+## Invalid bounds.
 ##
-## @table @asis
-## @item 204 (@w{@code{LPX_E_FAULT}})
-## Unable to start the search.
+## @item 5 (@w{@code{GLP_EFAIL}})
+## Solver failed.
 ##
-## @item 205 (@w{@code{LPX_E_OBJLL}})
+## @item 6 (@w{@code{GLP_EOBJLL}})
 ## Objective function lower limit reached.
 ##
-## @item 206 (@w{@code{LPX_E_OBJUL}})
+## @item 7 (@w{@code{GLP_EOBJUL}})
 ## Objective function upper limit reached.
 ##
-## @item 207 (@w{@code{LPX_E_ITLIM}})
+## @item 8 (@w{@code{GLP_EITLIM}})
 ## Iterations limit exhausted.
 ##
-## @item 208 (@w{@code{LPX_E_TMLIM}})
+## @item 9 (@w{@code{GLP_ETMLIM}})
 ## Time limit exhausted.
 ##
-## @item 209 (@w{@code{LPX_E_NOFEAS}})
-## No feasible solution.
+## @item 10 (@w{@code{GLP_ENOPFS}})
+## No primal feasible solution.
+##
+## @item 11 (@w{@code{GLP_ENODFS}})
+## No dual feasible solution.
+##
+## @item 12 (@w{@code{GLP_EROOT}})
+## Root LP optimum not provided.
 ##
-## @item 210 (@w{@code{LPX_E_INSTAB}})
+## @item 13 (@w{@code{GLP_ESTOP}})
+## Search terminated by application.
+##
+## @item 14 (@w{@code{GLP_EMIPGAP}})
+## Relative MIP gap tolerance reached.
+##
+## @item 15 (@w{@code{GLP_ENOFEAS}})
+## No primal/dual feasible solution.
+##
+## @item 16 (@w{@code{GLP_ENOCVG}})
+## No convergence.
+##
+## @item 17 (@w{@code{GLP_EINSTAB}})
 ## Numerical instability.
 ##
-## @item 211 (@w{@code{LPX_E_SING}})
-## Problems with basis matrix.
-##
-## @item 212 (@w{@code{LPX_E_NOCONV}})
-## No convergence (interior).
+## @item 18 (@w{@code{GLP_EDATA}})
+## Invalid data.
 ##
-## @item 213 (@w{@code{LPX_E_NOPFS}})
-## No primal feasible solution (LP presolver).
-##
-## @item 214 (@w{@code{LPX_E_NODFS}})
-## No dual feasible solution (LP presolver).
+## @item 19 (@w{@code{GLP_ERANGE}})
+## Result out of range.
 ## @end table
 ##
 ## @item extra
@@ -448,9 +429,28 @@
 ## @item time
 ## Time (in seconds) used for solving LP/MIP problem.
 ##
-## @item mem
-## Memory (in bytes) used for solving LP/MIP problem (this is not
-## available if the version of @sc{glpk} is 4.15 or later).
+## @item status
+## Status of the optimization.
+##
+## @table @asis
+## @item 1 (@w{@code{GLP_UNDEF}})
+## Solution status is undefined.
+##
+## @item 2 (@w{@code{GLP_FEAS}})
+## Solution is feasible.
+##
+## @item 3 (@w{@code{GLP_INFEAS}})
+## Solution is infeasible.
+##
+## @item 4 (@w{@code{GLP_NOFEAS}})
+## Problem has no feasible solution.
+##
+## @item 5 (@w{@code{GLP_OPT}})
+## Solution is optimal.
+##
+## @item 6 (@w{@code{GLP_UNBND}})
+## Problem has no unbounded solution.
+## @end table
 ## @end table
 ## @end table
 ##
@@ -481,7 +481,7 @@
 ## Author: Nicolo' Giorgetti <giorgetti@dii.unisi.it>
 ## Adapted-by: jwe
 
-function [xopt, fmin, status, extra] = glpk (c, A, b, lb, ub, ctype, vartype, sense, param)
+function [xopt, fmin, errnum, extra] = glpk (c, A, b, lb, ub, ctype, vartype, sense, param)
 
   ## If there is no input output the version and syntax
   if (nargin < 3 || nargin > 9)
@@ -608,7 +608,7 @@
     param = struct ();
   endif
 
-  [xopt, fmin, status, extra] = ...
+  [xopt, fmin, errnum, extra] = ...
     __glpk__ (c, A, b, lb, ub, ctype, vartype, sense, param);
 
 endfunction