David Avis          avis@cs.mcgill.ca     http://cgm.cs.mcgill.ca/~avis

  • What's New in Version 6.2
  • What's New in Version 6.1
  • What's New in Version 6.0
  • What's New in Version 5.0
  • What's New in Version 4.3
  • What's New in Version 4.2
  •  What's New in Version 4.1
  • What's New in Version 4.0

  • Introduction
  • lrslib: installation and usage
  • plrs: installation and usage
  • mplrs:installation and usage  (new in V6.0)

  • File formats
  • Basic options
  • Arithmetic packages 
  • Estimation    (revised in V6.0)

  •  Fourier Elimination
  • Linear Programming
  •  Nash Equilibria   (major revision in V6.1)
  • Volume and triangulation
  • Voronoi Diagrams
  • Extreme point enumeration and eliminating redundant inequalities
  • Linearities
  • Timing and interrupts
  • Error messages and troubleshooting
  • Hints and comments
  • Acknowledgements and References

    What's New in Version 6.2


  • mplrs:   new command-line options, improved performance on large-output instances, volume output added

  • What's New in Version 6.1

  • lrs: memory leaks fixed
  • mplrs: checkpointing after restart fixed
  • nash: new version of nash driver and new nash library, mostly by Terje Lensberg, memory leak and earlier bug in nash fixed

  • What's New in Version 6.0

  • mplrs: C wrapper prepared by Skip Jordan for lrs that allows for parallelization using the MPI library over a network of multi-core machines. It is derived from plrs but includes many modifications to ensure load balancing. It has delivered near linear speedups up to 1200 cores on runs by Kazuki Yoshizoe on Tsubame2.
  • Installation and usage are described here.
  • plrs is designed for shared memory machines and should improve on lrs even on dual core machines. mplrs is designed for machines and/or networks with at least 4 processors.
  •  Some additonal options were added to lrs in order to help with load balancing, in particular maxc

  • What's New in Version 5.0

  • plrs: C++ wrapper prepared by Gary Roumanis for lrs that allows for parallelization on multi-core machines and uses the Boost library. Installation and usage are described here and a technical description is described in this paper.
  •  Some changes were made to lrs in order to remove global variables and direct I/O, but its original functionality should be intact. Please report any missing functions.
  • The default arithmetic library has been changed to gmp, as this is distributed with most gcc installations. make all gives a gmp installation make allmp gives an installation with native lrs extended precision arithmetic
  • What's New in Version 4.3


  • Some corrections to the truncate option, which will produce all lex-positive bases representing a given vertex. See here for a description
  • What's New in Version 4.2c


  • 2nash is a 2-processor parallel version of nash                      (added 2009.2.18)
  • Bug fix for fourier: linearities disabled                                     (added 2009.2.18)
  • Triangulation of V-representation output via "volume" and "verbose" options
  • Estimator provides running time estimate
  • Bug fix related to memory allocation of linearity option.
  • Bug fix for nash, and inclusion of polytope version via setnash2.
  • incidence option now compatible with printcobasis n.

    What's New in Version 4.2a

  • Bug fixed relating to mis-scaled lp dual variables output when lponly set

    maxoutput n   Option limits output lines to n: either rays+vertices, or facets.

    What's New in Version 4.2

  • The main additions in Version 4.2 are two new drivers:


    There is a new bound option in lrs that truncates the reverse search tree when the objective function is less than the given bound for a maximization problem, or greater than the bound for minimization. When lrs is used to solve LPs, the dual variables will now be printed for the optimum solution.



    Other new basic options: dualperturb, project, printslack
    The nonnegative option should now work correctly  even if the origin is not a vertex of the polyhedron.


    What's New in Version 4.1

    Version 4.1 of lrslib contains no substantial revisions to the basic algorithms, but contains complete garbage collection and some simple drivers for solving vertex enumeration, convex hull and linear programming problems on generate polyhedra. The generation of input internally is simplified by the inclusion of some lpsolve like procedures to enable easy construction of the constraint matrix. Documentation for this is contained in http://cgm.cs.mcgill.ca/~avis/C/lrslib/lrslib.html

    A new feature of the basic package is the nonnegative option that speeds up the solution of problems with H-representation:

    b+Ax >= 0,        x>=0

    Note: The origin must be a vertex of the feasible region for Version 4.1.

    What's New in Version 4.0 (included in Version 4.1)

    Version 4.0 of lrs is distributed as the callable library lrslib. Two drivers are supplied, lrs.c and redund.c , that replace the corresponding programs in previous distributions. The new programs should be compatible with earlier ones, and correctly process existing input files. New in this version is that input files need not be of full dimension, and input equations or linearities may be specified. In addition, input files may contain redundant columns, which are closely connected to linearities. See also the section File Formats for more information.

    This distribution also comes with three  arithmetic lpackages: lrsmp, lrslong and lrsgmp . The library lrsmp is essentially the same extended precision package included in earlier distributions. lrslong is an implementation using long integer arithmetic, which is considerably faster, but does no overflow checking. lrsgmp is an interface to GNU MP which must be installed first. Binaries are produced for all arithmetic packages using the distributed make file.

    Most of the functions required for lrs and redund are included in the library lrslib. This is compiled with either of the arithmetic packages. One of the purposes of supplying a library is to allow simple customization of the reverse search procedure. The program lrs.c is the reverse search driver, and supplies each output vector to the user. It can be customized to prune the search according to various criteria, eg., if the value of some objective function falls below some value, or if an integer vector is produced etc.

    The installation procedure has been simplified, and the package is distributed as a compressed tar file. Included is a makefile for various configurations and some typical input files. Binaries are also included for some platforms.

    There are new options:  incidence , which lists all inequalities (resp. vertices/rays) which are incident with the current output vertex/ray (resp. facet). and truncate which prunes the search tree everytime a vertex different from the starting vertex is reached.


    Introduction

    A polyhedron can be described by a list of inequalities (H-representation) or as by a list of its vertices and extreme rays (V-representation).lrs is a C program that converts a H-representation of a polyhedron to its V-representation, and vice versa.  These problems are known respectively at the vertex enumeration and convex hull problems.
    Fukuda's FAQ page   contains a more detailed introduction to the problem, along with many useful tips for the new user.

    lrs is based on the reverse search algorithm of Avis and Fukuda(1992) , modified to use lexicographic pivoting  and implemented in rational arithmetic. See Avis(1998a) for a technical description, and Avis(1998b) for some computational experience. The input files are in Polyhedra format , developed by Fukuda and the author. The format is essentially self-dual, and the output file produced can be read in as an input file, with very minor modifications, to perform the reverse transformation. This format is compatible with that  used in Fukuda's cdd package, which performs the same transformations using a version of the double description methodcdd can also be used  in conjunction with lrs as a pre-processor for projections to subspaces, or as a post processor for computing the entire face lattice. Another program using the same file format is the primal-dual method pd, developed by Bremner, Fukuda and Marzetta .  It is essentially dual to lrs, and is very efficient for computing H-representations of simple polyhedra, and V-representations of simplicial polyhedra. It will compute the volume of a polytope given by an H-representation.

    Another program  based on the double description method is Christof and Loebel's porta , and a versatile tool for the algorithmic treatment of polytopes is Gawrilow and Joswig's  polymake package. Barber et al.'s qhull.  A package for volume computation called Vinci has been developped by Enge. A comprehensive general source of related infomation are Erickson's Computational Geometry Pages.

    In version 4.0, polyhedra handled by lrs  need not be full dimensional  and may contain input linearities and redundant columns .  lrs accepts either integer or rational input, and produces integer or rational output. All computations are done exactly using either extended precision arithmetic or fixed long integer arithmetic. In the latter case no overflow checking is performed, and the user is advised to check results using the extended precision version. Since it is a pivot based method, lrs can be very slow for degenerate inputs: i.e.. H-representations of non-simple polyhedra, and V-representations of non-simplicial polyhedra. On the other hand, it does not store the vertices/ rays or facets produced, so for very large problems it may be the only method that can solve the problem.  A discussion of various vertex enumeration/convex hull methods and the types of polyhedra that cause them to behave badly is contained in Avis, Bremner and Seidel( 1997).

    Additional functions of lrs include:

    A second program, redund, is supplied for removing redundancy from and H or V-representation. In the former case, this involves the removal of any inequalities that are not required to represent the polyhedron. In the latter case, this is  the problem of evaluating the extreme points and rays of a V-representation. These problems are normally  considerably easier than the transforamtions performed by lrs. In some cases, redundancy can greatly slow the processing time taken by lrs, and it is advisble to remove any redundancy from the input file using redund before applying lrs.

    These programs can be distributed freely under the GNU GENERAL PUBLIC LICENSE. Please read the file COPYING carefully before using.  Please inform the author of any interesting applications for which lrs/redund were helpful.


    lrslib installation and usage


    File formats

    Note for cdd users: lrs uses essentially the same file format as cdd. Files prepared for cdd should work with little or no modification. Note that  the V-representation corresponds to the "hull" option in cdd. Options specific to cdd can be left in the input files and will be ignored by lrs.  Note the input files for lrs are read in free format, after the line m n rational lrs will look for exactly m*n rationals or integers separated by white space (blank,  carriage return, tab etc.). lrs will not "drop" extra columns of input if n is less than the number of columns supplied.

    Basic options

    allbases bound  x                                    // Use with H-representation  - for lrs or nash //
    Either the maximize or minimize option should be selected. x is an integer or rational.
    For maximization (resp. minimization) the reverse search tree is truncated  whenever the current objective value is less (resp. more) than x .

    cache n debug  startingbasis endingbasis
    Print out cryptic but detailed trace, dictionaries etc. starting at #B=startingbasis and ending at #B=endingbasis. debug 0 0 gives a complete trace.
    digits n                // placed before the begin statement// dualperturb
    If lrs is executed with the maximize or minimize option, the reverse search tree is rooted at an optimum vertex for this function.
    If there are mulitiple optimum vertices, the output will often not be complete. This option gives a small perturbation to the objective to avoid this.
    A warning message is given if the starting dictionary is dual degenerate .

    estimates k geometric                 // H-representation  or voronoi option only //               This indicates ray  0   r0   r 1 ...   rn-1is incident with vertex  1   v0   v 1 ...   vn-1.  For the quadrant, the output is: incidence                                         
    This option automatically switches on printcobasis , so see below for a description of this option first.
    Can be used with printcobasis n. (Ver 4.2b)

    For input H-representation, indices of all input inequalities that contain the vertex/ray that is about to be output. For a simplicial face, there is no new output, since these indices are already listed. Otherwise, the additional tight inequalities are listed after a colon. Eg:
    V#1 R#0 B#1 h=0 facets  12 14 15 16 : 9 10 11 13 I#8 det= 8
     1  0  0  0  1

    The vertex 0 0 0 1 satisfies 8 input inequalities as equations, as indicated by I#8 : those with indices 12,14,15,16 are in the cobasis, and those with indices 9, 10, 11, 13 are in the basis. For a ray:
    V#1 R#5 B#1 h=0 facets  5 9* 10 11 12 13 : 2 3 4 I#8 det= 8
     0  1  1  0  0  1  1

    Here the ray 1  1  0  0  1  1 lies on 8 inequalities, with indices 5 10 11 12 13 in basis and 2 3 4 in cobasis. The starred index 9* indicates that the ray is terminated by the input inequality 9. This inequality is in the cobasis and defines the vertex from which the ray starts.

    For input V-representation, indices of all input vertices/rays that lie on the facet that is about to be output:
    F#5 B#3 h=2 vertices/rays  7 8* 11 13 15 : 1 3 5 9 I#8 det= 16
    1 -1  0  0  0

    The facet generated by inequality x1 <= 1 contains 8 input vertices, as indicated by I#8: those with indices 7,11,13,15 are in the cobasis, and those with indices 1 3 5 9 are in the basis.The starred index 8* indicates that this vertex  is also in the cobasis, but is not contained in the facet. It arises due to the lifting operation used with input V-representations.

    #incidence
    The same as printcobasis. Included for compatability with cdd.
    linearity  k  i1 i2 i ... ik
    The input contains k linearities in rows i1 i2 i ... ik of the input file are equations. See Linearities.
    maxdepth k maximize  a0 a1 ... an-1                                            // H-representation  only //
    minimize   a0 a1 ... an-1                                           // H-representation  only //
    If used with lrs the starting vertex maximizes (or minimizes) the function  a0 + a1 x 1 + ... + an-1 xn-1.
    The dualperturb option may be needed to avoid dual degeneracy.
    See Nash Equilibria and  Linear Programming
    maxcobases n         //from Version 6.0 //
           After n cobases have been generated lrs terminates and returns restart data for all unexplored roots of subtrees (except for leaves which are output). These subtrees are the unexplored siblings on the path back to the root of the reverse search tree. Used by mplrs to break up large subtrees into smaller pieces.

    maxoutput n
       
          
    Limits number of output lines produced (either vertices+rays or facets) to n

    mindepth k
    nonnegative                      // This option must come before the begin statement//
                                                                                                //H-representation only //
               Bug: Can only be used if the origin is a vertex of the polyhedron 
    For problems where the input is an H-representation of the form b+Ax>=0, x>=0 (ie. all variables non-negative, all constraints inequalities) it is not necessary to give the non-negative constraints explicitly if the nonnegative option is used. This option cannot be used for V-representations, or with the linearity option (in which case the linearities will be treated as inequalities). This option may be used with redund , but the implied nonnegativity constraints are not tested themselves for redundancy. To test everything it is necessary to enter the nonnegativity constraints explicitly in the input file. (In Ver 4.1, the origin must be a vertex).


    printcobasis  k                                  

    printslack                        // Use with H-representation //

    lrs prints a list of the indices of the input inequalities that are satisfied strictly for the current vertex, ie. corresponding slack variable is positive.
    If nonnegative is set, the list will also include indices n+i for each decision variable xi which is positive.


    project
    Used by driver fourier only.

    restart  V# R# B# depth {facet #s or vertex/ray #s}                startingcobasis i1 i2 i ... in-1 truncate                                            // H-representation only //      Some corrections made in v4.3
    The reverse search tree is truncated(pruned)  whenever a new vertex is encountered. Note: This does note necessarily produce the set of all vertices adjacent to the optimum vertex in the polyhedron, but just a subset of them. See here for a description of how to use this option.
    verbose
    Print slightly more detailed information about the run.
    volume                                              // V-representation  only // voronoi                                              // V-representation  only - place immediately after end statement //


    Arithmetic packages               

    The default arithmetic package in lrslib is GMP. This is usually available with gcc, but if not it needs to be installed separately from https://gmplib.org/
    In the lrslib makefile this library appears as:

      lrsgmp   An interface to GNU MP  which must be installed first.
     
    lrslib also has two build in packages that do not need separate installation, and can be build by % make allmp                

      lrsmp     Extended precision arithmetic 

     lrslong   Fixed length long integer arithmetic. No overflow checking
                       It is about 5-10 times faster than lrsmp and 4-5 times faster than lrsgmp.

      The standard make gives binaries  lrs/redund with lrsgmp, and lrs1/redund1 with lrslong
      Since no overflow checking is performed for lrs1/redund1 results should be confirmed by using one of the other arithmetic packages.

     


    Estimation

    maxdepth d
    estimates k            (New in V6.0) In order to try to get a better estimate, we implemented the additional option which makes the search go deeper than maxdepth when the estimates are large.

    subtreesize n     (default is MAXD)

                If an estimate obtained for a subtree at depth d exceeds n then the estimation is continued deeper into that subtree until an estimate less than or equal to n is obtained.

               The running time of lrs is proportional to the number of bases, so an estimate of the number of bases gives an easy way to estimate the running time for solving the complete problem by lrs:
               total running time = time for estimate * estimated number of bases / tree nodes evaluated.
    seed n lrs n < filename
     
      plrs

    C++ wrapper prepared by Gary Roumanis for lrs that allows for parallelization on multi-core machines and a technical description is described in this paper.
    Note: plrs performs only the standard lrs operation of converting a H-representation to V-representation or vice versa.

    Default installation:
    %make plrs                   (See below for more information)

    Usage is % plrs <infile> [ <outfile> ] [ -mt <max threads> ] [ -id <initial depth> ]

    -mt <max threads> specifies the number parallel threads calling lrs (default 12)

    -id <initial depth> specifies the initial depth of the RS tree to generate before parallelization (default 5)

    Note: Setting the two parameters correctly is very important ! Firstly it usually does not help to set -mt to a larger number than the number of cores on your machine.
    The parameter -id sets the initial depth for a single lrs run to populate a cobasis list for the parallelization. If -id is too small, the list will be too few to balance the load. If it is too large time is wasted on the initial non-parallel step. Furthermore, each cobasis on the list has to be reinitalized (ie. using the restart option in lrs) and this creates additional overhead. These issues are described in detail in Section 5 of the paper.

    Example: Input file mp5.ine is run with 8 threads after an initial descent to depth 3 in the tree.
    This produced 122 subtrees that were enumerated in parallel using 8 cores. The output file mp5.plrs is in the distribution.

    mai20%  plrs mp5.ine mp5.plrs -mt 8 -id 3
    *plrs:lrslib v.6.0 2015.7.13(lrsgmp.h)8 processes
    *Copyright (C) 1995,2015, David Avis   avis@cs.mcgill.ca
    *Input taken from mp5.ine
    *Output written to: mp5.plrs
    *Max depth of 3 to initialize starting cobasis list
    *Finished initializing cobasis list with 122 starting cobases
    *Starting 8 producer thread(s) and 1 consumer thread
    *Phase 1 time: 0 seconds
    *Totals: vertices=32 rays=0 bases=9041 integer-vertices=16
    *Elapsed time: 0 seconds
    :

    Installation instructions for plrs.

    1. Install boost library from http://www.boost.org/ into prefix/boost151 If you have root permission, prefix=/usr/include However you can install boost locally wherever you like. Instructions for installing the library are located here http://www.boost.org/doc/libs/1_51_0/more/getting_started/unix-variants.html. Look at section 5 for an easy install. ***Important: make a note of the path given at the end of the install process ****
    2. After the boost library is installed merge the contents of boost-atomic directory (included in lrslib-050) with the corresponding directories in your boost prefix/include/boost151 directory.
    3. Update the makefile to include the paths you recorded in step 1. "make plrs" will make plrs with the gmp library (assuming the gmp library is already installed on that machine). "make plrsmp" will make plrs with the standard lrsmp arithmetic library and plrs1 with the long integer library.
    4. Binaries using the native lrsmp and lrslong arithmetic libraries can be made by
    % make plrsmp                          (See arithmetic packgages for more details)

    mplrs

    C wrapper prepared by Skip Jordan for lrs that allows for parallelization using the MPI library over a network of multi-core machines. It is derived from plrs but includes many modifications to ensure load balancing. Near linear speedups up to 1200 cores have been observed. Generally plrs should be faster than mplrs when the number of cores are few, say 8 or less.
    Note: mplrs performs only the standard lrs operation of converting a H-representation to V-representation or vice versa.

    Default installation:
    % make mplrs                            (See below for more details)



    Default usage:    <number of processes> should be 4 or higher

    % mpirun -np <number of processes>  mplrs <infile> [ <outfile> ]     (See below for all mplrs options)

    Note: Remove all options after the end statement of the input file!

    Example: Input file mp5.ine is run with 8 processors. The output file mp5.mplrs is in the distribution.
    This produced 378 subtrees that were enumerated in parallel using 6 producer cores,  1 core controlling the run and 1 core collecting the output.

    mai20% mpirun -np 8 mplrs mp5.ine mp5.mplrs
    *mplrs:lrslib v.6.0 2015.7.13(lrsgmp.h)8 processes
    *Copyright (C) 1995,2015, David Avis   avis@cs.mcgill.ca
    *Input taken from mp5.ine
    *Output written to: mp5.mplrs
    *Starting depth of 2 maxcobases=50 maxdepth=0 lmin=3 lmax=3 scale=100
    *Phase 1 time: 0 seconds.
    *Total number of jobs: 378, L became empty 4 times
    *Totals: vertices=32 rays=0 bases=9041 integer-vertices=16
    *Elapsed time: 1 seconds.
    2.285u 0.137s 0:01.86 129.5%    0+0k 0+9976io 36pf+0w





    Installation instructions for mplrs.


    1. Install MPI library on all machines that will be used for a single run. We have tested mplrs with both openmpi and mpich2. Note that if you install MPI from binary packages, you will need both the library and development packages if available. Configuring MPI to run on a cluster of machines is beyond the scope of this document, however mplrs should work well on any properly-configured cluster or supercomputer.
    2. Make any path changes to the beginning of makefile as necessary.
    3. Run %make mplrs
    Note 1: exactly the same version of mplrs must be used on each machine in the cluster, or unexpected results may occur!
    Note 2: mplrs creates temporary files on each machine in the cluster. The default is to use /tmp. If this is not writeable the -temp option must be used, see below.
    4. A binary using the lrslong fixed integer arithmetic library can be made by running
    % make mplrs1                           (See arithmetic packgages for more details on lrslong)

    Full set of options for mplrs
    mplrs
    has many options, but just using the default settings should provide good performance in general.
    The parameters -hist and -freq give interesting information about the degree of parallelization.

     mpirun -np <number of processes> mplrs <infile> <outfile> -id <initial depth> -maxc <maxcobases> -maxd <depth> -lmin <int> -lmax <int> -scale <int> -hist <file> -temp <prefix> -freq <file> -stop <stopfile> -checkp <checkpoint file> -restart <checkpoint file> -time <seconds>
                                                              
     -np <number of processes>     required, should be at least 4
                                                         default
     -id <initial depth>                      2                      the depth of the original tree search to populate the job queue L
    -maxc <maxcobases>                50 (*scale)     a producer stops and returns all subtrees that are not leaves to L after generating  maxc nodes   
    -maxd <depth>                            0                       a producer returns all subtrees that are not leaves at depth maxd. Zero if not used
    -lmin <int>                                    3                       if job queue |L|<np*lmin then the maxd parameter is set for all producers         
    -lmax <int>                                   lmin                 if job queue |L|>np*lmax then then maxc is replaced by maxc*scale
    -scale <int>                                   100                   used by lmax
    -hist <file>                                                               store parallelization data in <file> for use by gnuplot, see below
    -temp <prefix>                          /tmp/                  store a temporary file for each process. Should be specified if /tmp not writeable
    .
                                                                                       Using " -temp  ./ " will write temporary files to the current directory
    -freq <file>                                                              store frequency data in <file> for use by gnuplot, see below
    -stop <stopfile>                                                     terminate mplrs if a file with name <stopfile> is created in the current directory
    -checkp <checkpoint file>                                  if mplrs is terminated by -stop or -time then it can be restarted using this <checkpoint file> and -restart
    -restart <checkpoint file>                                  restart mplrs using previously created <checkpoint file>. If used with -checkp file names should be different!
    -time <seconds>                                                   terminate mplrs after <seconds> of elapsed time

    New in V6.2:
    -countonly                                                             don't output vertices/rays/facets, just count them
    -maxbuf <n>                              500                   controls maximum size of worker output buffers
                                                                                     smaller values increase "streaminess" of the output while larger values send larger blocks of output
    -stopafter <n>                                                     exit after approximately <n> cobases have been computed (no guarantee about how many vertices/rays/facets computed)

    Checkpoints, restarts and abnormal termination

    mplrs has inherited restart capability from lrs, which was described in Timing and interrupts
    Details of the layout of the checkpoint file can be found here.
    Checkpointing can be very useful for long runs that use significant computational resources. If resources are needed elsewhere, or additional resources become available, it is possible to stop mplrs and use the checkpoint file to restart. The mpirun command can be modified for this restart to reflect the new resource availability. The normal way to do this by using the -checkp and either the -stop or -time options. In addition if mplrs receives a SIGTERM or SIGHUP signal, it checkpoints and terminates.  This allows users to checkpoint runs that were not started with -stop <file> or -time <seconds>.  These signals can be sent using the kill utility:

    % kill <pid>
    where <pid> is the process ID of an mplrs instance to be terminated.
    Sending the signal to any one of the mplrs instances is sufficient.

    The checkpoint is produced in the file specified by the -checkp
    argument.  If it was not specified, the checkpoint is produced in the output.
    In this case, the checkpoint file is contained between the lines
     *Checkpoint file follows this line
    and
     *Checkpoint finished above this line .
    These two lines are not part of the checkpoint file.

    To terminate the run immediately, without a checkpoint, use
    % kill -9 <pid>

    Visualization of parallelization
    A good realtime view of an mplrs run can be obtained by using the -hist and -freq parameters and gnuplot using the plotL.gp and plotD.gp files.
    These files can be used if the parameters are set as  -hist hist -freq freq and can be obtained while mplrs is running as follows:
    %gnuplot plotL.gp
    %gv plotL.ps
    %gnuplot plotD.gp
    %gv plotD.ps
    If other files are used for -hist and -freq then plotL.gp and plotD.gp should be edited accordingly.
    The first plot has 3 graphs showing the number of processors working, the size of the job queue L, and message requests pending, all versus elapsed time on the x-axis.
    The second histogram shows the size of the subtrees explored by producers. The root of the subtree is not counted in this size. Note that a producer stops after maxc nodes are explored, but in backtracking some additional leaves may be discovered. So the size of the largest subtree maybe slightly larger than maxc.

    Details of -hist and -freq files

    mplrs periodically adds a line to the histogram file, for example
       54.118141 94 279 94 0 0 373

    This line contains the following information:
       Time since execution began in seconds (54.118141 here)
       Number of busy workers (94 here)
       Current size of job queue (279 here)
       Number of workers that may return unfinished jobs (94 here)
       Unused (0 here)
       Unused (0 here)
       Total number of jobs that have been in the job queue (373 here)
    The second and fourth entries are similar and can differ due to
    latencies.

    The frequency data file contains one integer per line, with each integer
    corresponding to the size of a subtree explored by a worker.



    Fourier Elimination (Due to several bug reports, this program should not be used. Caveat emptor!)

    Tallman Nkgau has contributed a driver fourier that computes the projection of a polyhedron given by its H-representation onto a selected set of coordinates. The program can be compiled by

    make fourier

    and is run by the command

    fourier file.ine [fileout]

    file.ine
    is a standard H-representation that does not contain linearities. If you have linearities, please replace each by two inequalities.
    Following the end statement, insert the option:

    project t  a1 ... at

    The output will be the polyhedron projected onto the coordinates a1 ... a , each of which is a unique number betwen 1 and n.


    Nash Equilibria                                              Major revision in V6.1             

    See:    D. Avis, G. Rosenberg, R. Savani, B. von Stengel, "Enumeration of Nash Equilibria for Two-Player Games", Economic Theory 42(2009) 9-37  pdf

    lrsnash computes all Nash equilibria (NE) for a two person noncooperative game are computed using two interleaved reverse search vertex enumeration steps.
    The input for the problem are two m by n matrices A,B of integers or rationals. The first player is the row player, the second is the column player.
    If row i and column j are played, player 1 receives Ai,j and player 2 receives Bi,j.

    Version 6.1 contains lrsnash.c and lrsnashlib.c replacing nash.c with a library version and simpler interface that does not require setupnash. Big thanks to Terje Lensberg for this.
    nashdemo.c is a very basic template for setting up games and calling the library function lrs_solve_nash(game *g).

    To compile:

    % make lrsnash

    you get binaries lrsnash, nashdemo and 2nash.

    % nashdemo                     just runs the demo, no parameters

    % lrsnash game                 finds the equlibrium for file game

     [ % lrsnash game1 game2      or    % 2nash game1 game2               finds the equilibrium for game1 and game2 in legacy nash format which is  described below ]

    The file game is in the format :
     m n
    matrix A
    matrix B

    eg. the file game is for a game with m=3 n=2:
    3 2

    0 6
    2 5
    3 3

    1 0
    0 2
    4 3



    % lrsnash game

    *Copyright (C) 1995,2015, David Avis   avis@cs.mcgill.ca
    *lrsnash:lrslib v.6.1 2015.10.27(lrsgmp.h gmp v.6.0)
    2  1/3  2/3  4
    1  2/3  1/3  0  2/3

    2  2/3  1/3  3
    1  0  1/3  2/3  8/3

    2  1  0  3
    1  0  0  1  4

    *Number of equilibria found: 3
    *Player 1: vertices=5 bases=5 pivots=8
    *Player 2: vertices=3 bases=1 pivots=8

    *lrsnash:lrslib v.6.1 2015.10.27(32bit,lrsgmp.h)

    ------------------------------------------------------------------
    Output interpretation:
    Each row beginning 1 is a strategy for the row player yielding a NE with each row beginning 2 listed immediately above it.
    The payoff for player 2 is the last number on the line beginning 1, and vice versa.
    Eg: first two lines of output: player 1 uses row probabilities 2/3 2/3 0 resulting in a payoff of 2/3 to player 2.
    Player 2 uses column probabilities 1/3 2/3 yielding a payoff of 4 to player 1.


    Legacy format

    % setnash game game1 game2    

    produces two H-representations, game1 and game2, one for each player.

    To get the equilibrium, run

    % lrsnash game1 game2

    If you have two or more cpus available run 2nash instead of lrsnash as the order of the input games is immaterial.
    It runs in parallel with the games in each order. (If you use lrsnash, the program usually runs faster if m is <= n , see below.)

    % 2nash game1 game2 [outputfile]

    If no output file is specified the output is placed in a file names out.

    If both matrices are nonnegative and have no zero columns, you may instead use setnash2:

    % setnash2 game game1 game2

    Now the polyhedra produced are polytopes.

    The output  of lrsnash in this case is a list of unscaled probability vectors x and y.

    ----------------------------------------------------------
    % lrsnash game1 game2
    Processing legacy input files. Alternatively, you may skip
    setupnash and pass its input file to this program.

    *nash:lrslib v.6.1 2015.10.27(lrsgmp.h gmp v.6.0)
    *Copyright (C) 1995,2015, David Avis   avis@cs.mcgill.ca
    *Input taken from file game1
    *Second input taken from file game2

    ***** 4 3 rational
    2  1/12  1/6
    1  1  1/2  0

    2  2/9  1/9
    1  0  1/8  1/4

    2  1/3  0
    1  0  0  1/4


    *Number of equilibria found: 3
    *Player 1: vertices=6 bases=6 pivots=8
    *Player 2: vertices=4 bases=1 pivots=8
    *nash:lrslib v.6.1 2015.10.27(32bit,lrsgmp.h)

    -----------------------------------------------------------


    To normalize, divide each vector by v = 1^T x and u=1^T y.
    u and v are the payoffs to players 1 and 2 respectively.





    If m is greater than n then lrsnash usually runs faster by transposing the players. This is achieved by running:

    % lrsnash game2 game1

    If you wish to construct the game1 and game2 files by hand, they are fragile and should be done exactly as follows:

       For player 1: eg. game1
       One linearity in the last row
       Identity matrix with additional final column 0
       Transpose of payoff matrix for player 2 with final column 1
       Last row is prob sum to one

       For player 2: eg. game2
       One linearity in the last row
       Payoff matrix for player 1 with final column 1
       Identity matrix with additional final column 0
       Last row is prob sum to one

    Corresponding to file game above we get

    *game: player 1
    H-representation
    linearity 1 6
    begin
    6 5 rational
    0 1 0 0 0
    0 0 1 0 0
    0 0 0 1 0
    0 -1 0 -4  1
    0 0 -2 -3  1
    -1 1 1 1 0
    end

    *game: player 2
    H-representation
    linearity 1 6
    begin
    6 4 rational
    0 0 -6  1
    0 -2 -5  1
    0 -3 -3  1
    0 1 0 0
    0 0 1 0
    -1 1 1 0
    end





      Linear Programming             

    lponly

                 and one of the options maximize or minize:

    maximize a0 a1 ... an-1                                           // H-representation  only //

    minimize a0 a1 ... an-1                                             // H-representation  only //

    To print the dictionary at a few key points also include the option:

    verbose

    New in V4.2. Dual variables are now printed at termination. If the linearity option is used, only a partial list of dual variables will be given.
                           Dual variable yi refers to inequality number i in the input.


    Volume and triangulation

    lrs can be used to compute the volume of a full dimensional polytope given as a V-representation. This follows from the fact that lex-postive bases form a triangulation of the facets, and that a V-representation is always lifted. See "Theoretical Description" on lrs home page for some remarks on this. The option

    volume                                                                             // V-representation only //

    will cause the volume to be computed. For input cube.ext, the output is:
    *Volume=8

    The triangulation can be output by adding also the option verbose.
    This would give the output:

    F#0 B#1 h=0 vertices/rays  4 6 7 8 I#8 det= 8
     1  1  0  0
     1  0  1  0
     1  0  0  1
    F#3 B#2 h=1 vertices/rays  4 5 6 7 I#8 det= 8
    F#3 B#3 h=2 vertices/rays  3 4 5 7 I#8 det= 8
     1 -1  0  0
    F#4 B#4 h=3 vertices/rays  2 3 4 5 I#8 det= 8
     1  0  0 -1
    F#5 B#5 h=4 vertices/rays  1 2 3 5 I#8 det= 8
    F#5 B#6 h=2 vertices/rays  2 4 5 6 I#8 det= 8
     1  0 -1  0
    end
    *Sum of det(B)= 48
    *Volume= 8

    Each of the 6 bases corresponds to a simplex.
    The first simplex is composed of vertices 4 6 7 8, second simplex is 4 5 6 7, etc.

    If the volume option is applied to an H-representation, the results are not predictable. If the option is applied to a V-representation of  a polytope that is not full dimensional, the volume of a projected polytope is computed. The projection used is to the lexicographically smallest coordinate subspace, see Avis, Fukuda, Picozzi (2002)

    For polytopes given by a H-representation, it will first be necessary to compute the V-representation.


    Voronoi diagrams

    lrs can be used the compute the V-vertices of a Voronoi diagram of a set of data points in n-1 dimensional space. To do this we use a standard lifting procedure (see, e.g., Edelsbrunner, "Algorithms in Combinatorial Geometry," pp 296-297) . Each point is mapped to a half space tangent to the parabaloid in n dimensions, by the mapping:

    p , p , ...., p n-1     ->    (p1 +   p22 +  ...   +  pn-1 ) - 2 p1  x - 2 p x2 - .... - 2  p n-1 xn -1  + x n>= 0

    lrs is applied to the H-representation so created.  This transformation is performed automatically for a V-representation if the

    voronoi          // V-representation only - place immediately after end statement //

    option is specified.
    Note: The input file must consist entirely of data points (no rays), i.e.. there must be a one in column one of each line. The volume option should not be used, since the volume reported will not be the volume of the original V-representation.
    The output will consist of the Voronoi vertices (columns beginning with a one) and Voronoi rays (columns beginning with zero) for the Voronoi diagram defined on the data points.  If the printcobasis option is given, the n "data points" indices produced will tell which set of input data points corresponds to the given Voronoi vertex or ray. In case of degeneracies, a given Voronoi vertex may be generated by more than n of the input data points. In this case, use of the allbases option will cause all  sets of n input data points corresponding to a Voronoi vertex to be printed. For Voronoi rays, the immediately preceding  cobasis is the cobasis of the the Voronoi vertex from which the ray emanates.  The index followed by a * is the data point to drop in order to generate the ray. If the geometric option is given the correspondence between Voronoi rays and Voronoi vertices will be produced automatically.

    Example: Compute the Voronoi diagram of the planar point set (0,0), (2,1), (1,2), (0,4), (4,0), (4,4) (2,-4).
    vor7-3
    *6 Voronoi vertices and 5 rays
    *7 input data points
    V-representation
    begin
    7 3 integer
    1 0 0
    1 2 1
    1 1 2
    1 0 4
    1 4 0
    1 4 4
    1 2 -4
    end
    voronoi
    printcobasis
    allbases
    geometric

    The output produced is
     

    V-representation
    begin
    ***** 3 rational
    V#1 R#0 B#1 h=0 data points  1 5 7 det=64
    1 2 -3/2
    V#1 R#1 B#1 h=0 data points  1 5* 7 det=64
    0 -2 -1  * 1 2 -3/2
    V#1 R#2 B#1 h=0 data points  1* 5 7 det=64
    0 2 -1  * 1 2 -3/2
    V#1 R#2 B#2 h=1 data points  1 2 5 det=16
    1 2 -3/2
    V#2 R#2 B#3 h=2 data points  1 2 3 det=12
    1 5/6 5/6
    V#3 R#2 B#4 h=3 data points  1 3 4 det=16
    1 -3/2 2
    V#3 R#3 B#4 h=3 data points  1 3* 4 det=16
    0 -1 0  * 1 -3/2 2
    V#4 R#3 B#5 h=2 data points  2 5 6 det=32
    1 15/4 2
    V#4 R#4 B#5 h=2 data points  2* 5 6 det=32
    0 1 0  * 1 15/4 2
    V#5 R#4 B#6 h=3 data points  2 3 6 det=20
    1 27/10 27/10
    V#6 R#4 B#7 h=4 data points  3 4 6 det=32
    1 2 15/4
    V#6 R#5 B#7 h=4 data points  3* 4 6 det=32
    0 0 1  * 1 2 15/4
    end

    The output contains 6 Voronoi vertices :
    (2, -3/2), (5/6,5/6),(-3/2,2),(15/4,2), (27/10,27/10), (2,15/4).
    The Voronoi vertex (2,-3/2) appears twice in the output with data point indices 1 5 7 and 1 2 5. This means that it is degenerate and is defined by the set of 4 input data point in positions 1,2,5,7 in the input file. I.e.. it is the centre of an empty  circle through the four input data points (0,0), (2,1), 4,0), (2,-4).  The other Voronoi vertices appear once each and are defined respectively by the data points with indices (i.e..  position in the input file)  1 2 3,  1 3 4,  2 5 6,  2 3 6 and 3 4 6. The  Voronoi diagram has 5 rays
    (2, -3/2) + (-2t,-t),    (2,-3/2)+(2t,-t),    (-3/2,2)+(-t,0),    (15/4,2)+(t,0),    (2,15/4)+(0,t)

    For example, the first ray in the output appears:
    V#1 R#1 B#1 h=0 data points  1 5* 7 det=64
    0 -2 -1  * 1 2 -3/2

     This means that the ray (-2t,-t) emanates from the vertex defined by data points 1 5 7, namely (2, -3/2). The asterisk on index 5 indicates that the ray is defined by the data points with indices 5 and 7, namely (0,0) and (2,-4).
     
     


      Extreme point enumeration and eliminating redundant Inequalities

    A convex hull problem that occurs frequently is to enumerate the extreme points (vertices) of a given set of input points. This problem is in fact much simpler than the problem of finding the facets of the given input point set. It can be solved by linear programming.  The dual problem is to remove redundant inequalities from an H-representation. An input  inequality is redundant if it can be deleted without changing the polyhedron. The program redund also solves these problems. They can be solved by cdd using the vertex_listing  and facet_listing  options.

    To remove input points that are not vertices from a V-representation or redundant inequalities from an H-representation use the command:

    The resulting file can be used directly with lrs, or even piped into lrs. In fact, lrs works best if the input is non-redundant, see the section Redundancy vs Degeneracy.

    Note: Versions of redund prior to this release failed to remove redundancy from the starting basis.


    Linearities   


    linearity  k  i1 i2 i ... ik

    The input file contains k linearities. If the input is a H-representation, the rows i1 i2 i ... ik of the input file are equations. For a V-representation, the rows with these indices should begin with zero in column one, and will be interpreted as lines rather than rays.  Linearities defined on the input vertices of a V-representation are not defined, but the program will accept them and produce some output. Each of the indice ik must be a distinct number between 1 and m. With an  H-representation, linearities are useful for enumeration of vertices on a facet or lower dimensional subspace. For example the file:

    cube_ridge
    *cube of side 2 centred at the origin
    H-representation
    linearity 2  1 5
    begin
    6 4 rational
    1 1 0 0
    1 0 1 0
    1 0 0 1
    1 -1 0 0
    1 0 -1 0
    1 0 0 -1
    end

    causes vertices to be enumerated on the ridge which is the intersection of the two facets

    x1 = -1   and   x2 = 1

    so the output is the pair of vertices

    cube_ridge
    *Input linearity in row(s) 1 5
    V-representation
    begin
    2  4  rational
     1 -1  1  1
     1 -1  1 -1
    end

    Specifying linearities in this way will often produce redundancy , especially if the dimension of the problem is reduced considerably. As a preprocessing step, it is useful to apply to remove any redundancy by redund. In the case of the above problem the output produced by redund is:

    cube
    *Input linearity in row(s) 1 5
    *row 2 was redundant and removed
    *row 4 was redundant and removed
    H-representation
    linearity 2 1 2
    begin
    4 4 rational
     1  1  0  0
     1  0 -1  0
     1  0  0  1
     1  0  0 -1

    and two redundant halfspaces were removed.

    Redundant columns are closely related to linearities. If we examine the V-representation of cube_ridge above we can see that it is just a line segment in 3 dimensional space. Further,  columns 2 and 3 are multiples of column 1. If lrs is applied to this file, the column redundancies give rise to two linearities, so the output will appear as the H-representation given above: geometrically the intersection of two planes (the linearities) with two half-planes (defining the endpoints of the line segment).

    In general, the representation of the linearity space is not unique, however the one produced by lrs should be the same as that produced by cdd.



    Timing and interrupts

    lrs handles certain signals unless it is compiled with the -DOMIT_SIGNALS option. It is possible to interrupt lrs and get the latest cobasis, which can be used for restarting the program (useful if the machine is going down!)
            signal                     operation
            USR1                         print current cobasis and continue
            TERM                       print current cobasis and terminate
            INT (ctrl-C)                      ditto
            HUP                                      ditto
    lrs also provides timing information, unless compiled with the option -DOMIT_TIMES.

    Error messages and troubleshooting

    The most common error occurs from an incorrect input file specification, please check the section File Formats carefully. In particular, lrs does not check the type or number of input coefficients specified.  After the line
    m n rational
    you must specify exactly m*n rational or integer coefficients. They are read  in free format , but normally each input facet or vertex/ray is begun on a new line.  See note for cdd users.

    The following error messages are produced by lrs . They are  arranged in alphabetic order.

    Cannot find linearity in the basis

    The linearity option was specified but a basis cannot be created. Check the linearity indices are all less than n-1 and are disitinct.

    Data type must be integer of rational

    Digits must be at most 2295  Change MAX_DIGITS and recompile     (This message does not appear if the default gmp arithmetic package is used) Invalid input: check you have entered enough data!
    Usually means that end of file was reached before enough input data was read.

    Invalid Co-basis - does not have correct rank
    Maximize/minimize only valid for H-representation No begin line No data in file No feasible solution Starting cobasis indices must be distinct and in range 1 .. m Trying to restart from infeasible dictionary


    Hints and comments

    H- vs V- representation

     lrs is programmed to manipulate H-representations directly. A file presented as a V-representation is processed by lifting it to a cone in one higher dimension, which is treated internally as a H-representation. If the input file is a polytope which contains the origin, then the user has two options. Submit it as a V-representation and have it processed as just described, or submit it as a H-representation, and interpret the output as a list of facet inequalities rather than "vertices". Since this will not be lifted, it will be processed in a different way by lrs. Sometimes a degenerate V-representation may run more quickly as a H-representation, and sometimes more slowly. To decide which representation to use for a large problem, the user can run the estimates option and choose the representation with fewest estimated bases.
     

    Redundancy vs Degeneracy

    For an H-representation, an input is redundant if some inequality can be deleted without changing the polyhedron. It is degenerate if (in d dimensions) at least one vertex lies on d+1 or more facets.  Similarly in a V-representation an input is redundant if some input point is not a vertex of the convex hull.  It is degenerate if some facet contains d+1 or more input points. The options   printcobasis and incidence give degeneracy information. Degeneracy causes pivot  or triangulation based methods such as lrs to  run slowly. Redundancy is one cause of degeneracy, but it can be avoided by pre-processing the input files. See section Extreme Point Enumeration and Redundant Inequalities for instructions on how to do this. This pre-processing is unnecessary if it is known that the input is non-redundant.

    Even with redundant input removed a polyhedron may be highly degenerate. In directory ine/metric there are many highly degenerate combinatorial polytopes. These are difficult problems for all vertex enumeration/convex hull programs that use pivoting, such as lrs.  For example, the file ccc7.ine is a cone with 63 facets in 21 dimensions. It has 38,780 extreme rays, but computing these required the evaluation of 247,271,659 bases!

    Memory considerations

    The strong point of lrs is that it does not save the output produced, so in theory it cannot run out of memory.  With cache size one all memory is allocated at the beginning, so if lrs starts running it will not run out of memory. It is possible however that the number of digits required to do the calculations exceeds the amount specified on the digits option, or the default. In practice, this problem will also arise early in the computation. In any case, a message is printed and the calculation can be restarted. In order to improve performance, some dictionaries should be cached. The default of 10 can be overridden by the cacheoption. If the dictionary is in the cache it does not need to be recomputed when backtracking, reducing  processing time by about 40%. Since the cache is allocated dynamically, a cache size that is too large can potentially use up large ammounts of machine memory.

    Geometric Rays

    A minimum V-representation of a polyhedron is a minimum set of vertices and rays such that each point in the polyhedron can be expressed as a convex combination of vertices plus a non-negative combination of rays. For the cube, if we delete the inequality
    x3 <= 1, i.e.. the line 1 0 0 -1 from file cube.ine, we get the output:
    V-representation
    ***** 4 rational
    1 1 1 -1
    0 0 0 1
    1 -1 1 -1
    1 1 -1 -1
    1 -1 -1 -1
    end

    indicating the polyhedron is the convex combination of 4 vertices and 1 ray. With the geometric option, we get the output:
    V-representation
    begin
    ***** 4 rational
    1 1 1 -1
    0 0 0 1  * 1 1 1 -1
    1 -1 1 -1
    0 0 0 1  * 1 -1 1 -1
    1 1 -1 -1
    0 0 0 1  * 1 1 -1 -1
    1 -1 -1 -1
    0 0 0 1  * 1 -1 -1 -1
    end

    This indicates that geometrically, the polyhedron has 4 parallel extreme rays (0,0,t) , one incident to each vertex. With the geometric option, all rays will be printed. Without the option, lrs tries to print each ray once, but in some cases duplicates will remain, see  subsection Output Duplication.

    Output Duplication

    For degenerate inputs, pivot based methods for vertex/ray enumeration such as lrs may generate the same output ray many times. An output is only printed when it occurs with a lexicographically minimum basis. This removes all duplicate vertices, but rays may still be output more than once. This is due to the fact that duplicate geometric rays cannot always be detected without storing the output. Since V-representations are automatically lifted to a higher dimension, this will not happen for facet enumeration. Unless the allbases option is specified, lrs makes checks in order to remove duplicates.   A warning message is produced when duplicates may occur in the output. They can be removed using the program buffer.c. Two important types of input never produce duplicate output: polytopes (i.e. bounded polyhedra) and cones (i.e. polyhedra where the origin is the only vertex).


    Acknowledgements and References

    I would like to thank many people for helping with this implementation project. Komei Fukuda encouraged me from the start, collaborated in designing the file formats, and provided many suggestions for improving the code. Debugging would have been almost impossible without the use of his program cdd as a benchmark. David Bremner implemented memory allocation, cacheing and signals. Ambros Marzetta demonstrated the importance of cacheing and lrslong is based on his earlier implementation of this as prs_single.  Jerry Quinn coded the integer divide routine. Bug reports were provided by many users, for which I thank them. In particular Gerardo Garbulsky's extensive use of earlier versions suggested many refinements and Andreas Enge helped debug the volume computation. Tallman Nkgau contributed fourier.

    D. Avis, lrs: A Revised Implementation of the Reverse Search Vertex Enumeration Algorithm, http://cgm.cs.mcgill.ca/~avis/doc/avis/Av98a.ps
       In: Polytopes - Combinatorics and Computation, Ed. G. Kalai and G. Ziegler, Birkhauser-Verlag (2000) 177-198.

    D. Avis, "Computational Experience with the Reverse Search Vertex Enumeration Algorithm," Optimization Methods and Software, (1998 (to appear)). http://cgm.cs.mcgill.ca/~avis/doc/avis/Av98b.ps

    D. Avis, D. Bremner, and R. Seidel, "How Good are Convex Hull Algorithms?," Computational Geometry: Theory and Applications, Vol 7,pp.265-301(1997). http://cgm.cs.mcgill.ca/~avis/doc/avis/ABS96a.ps

    D. Avis and L. Devroye, "Estimating the Number of Vertices of a Polyhedron," pp. 179-190 in Snapshots of Computational and Discrete Geometry, ed. D. Avis and P. Bose, School of Computer Science, McGill University (1994). http://cgm.cs.mcgill.ca/~avis/doc/avis/AD94a.ps
      In: Information Processing Letters, (2000) V. 73, pp. 137-143.

    D. Avis and K. Fukuda, "A Pivoting Algorithm for Convex Hulls and Vertex Enumeration of Arrangements and Polyhedra," Discrete and Computational Geometry, Vol. 8, pp. 295-313 (1992).  http://cgm.cs.mcgill.ca/~avis/doc/avis/AF92b.ps

    D. Avis, K. Fukuda and S. Picozzi, "On Canonical Representations of Convex Polyhedra", Mathematical Software,  ICMS 2002, Ed. A. Cohen, X-S Gao, N. Takayama, World Scientific, pp.350-360 (2002)   http://cgm.cs.mcgill.ca/~avis/doc/avis/AFP02a.ps

    D. Avis, G. Rosenberg, R. Savani, B. von Stengel, "Enumeration of Nash Equilibria for Two-Player Games", Economic Theory 42(2009) 9-37  pdf

    D. Bremner, K. Fukuda and A. Marzetta, Primal-Dual Methods for Vertex and Facet Enumeration, 13th ACM  Symposium on Computational Geometry SCG 1997, 49-56.   http://www.cs.unb.ca/profs/bremner/pd/