cmu 1.7 → 1.8

raw patch · 3 files changed

+65/−57 lines, 3 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• ChangeLog +13/−0
• cmu.cabal +1/−1
• src/Algebra/CommutativeMonoid/LinDiophEq.hs +51/−56

@@ -1,3 +1,16 @@+2012-05-12  John D. Ramsdell  <ramsdell@mitre.org>++	* cmu.cabal (Version): Released as version 1.8.++	* src/Algebra/CommutativeMonoid/LinDiophEq.hs: Used the method in+	the paper by Contejean and Devie for solving an inhomogeneous+	equation using the algorithm for the homogeneous equation solver.+	The previous version buggy.++2012-05-02  John D. Ramsdell  <ramsdell@mitre.org>++	* cmu.cabal (Version): Released as version 1.7.+ 2012-04-27  John D. Ramsdell  <ramsdell@mitre.org>          * src/Algebra/CommutativeMonoid/LinDiophEq.hs: Corrected spelling

@@ -1,5 +1,5 @@ Name:			cmu-Version:		1.7+Version:		1.8 Maintainer:		ramsdell@mitre.org Cabal-Version:		>= 1.2 License:		GPL

@@ -25,63 +25,39 @@ -- The solver uses the algorithm of Contejean and Devie as specified -- in \"An Efficient Incremental Algorithm for Solving Systems of -- Linear Diophantine Equations\", Information and Computation--- Vol. 113, pp. 143-174, 1994 after a modification explained below.+-- Vol. 113, pp. 143-174, 1994. -- -- The algorithm for systems of homogeneous linear Diophantine -- equations follows.  Let e[k] be the kth basis vector for 1 <= k <= -- n.  To find the minimal, non-negative solutions M to the system of--- equations sum(i=1,n,a[i]*v[i]) = 0, the modified algorithm of--- Contejean and Devie is:+-- equations sum(i=1,n,a[i]*v[i]) = 0, the algorithm of Contejean and+-- Devie is: -- --  1. [init] A := {e[k] | 1 <= k <= n}; M := {} -- --  2. [new minimal results] M := M + {a in A | a is a solution} -----  3. [breadth-first search] A := {a + e[k] | a in A, 1 <= k <= n,--- \<sum(i=1,n,a[i]*v[i]),v[k]> \< 0}------  4. [unnecessary branches] A := {a in A | all m in M : some+--  3. [unnecessary branches] A := {a in A | all m in M : some --     1 <= k <= n : m[k] < a[k]} -----  5. [test] If A = {}, stop, else go to 2.+--  4. [breadth-first search] A := {a + e[k] | a in A, 1 <= k <= n,+-- \<sum(i=1,n,a[i]*v[i]),v[k]> \< 0} ----- The original algorithm reversed steps 3 and 4.+--  5. [test] If A = {}, stop, else go to 2. -- -- This module provides a solver for a single linear Diophantine -- equation a*v = b, where a and v are vectors, not matrices.--- Conceptually, it uses the homogeneous solver after appending -b as--- the last element of v and by appending 1 to a at each step in the--- computation.  The extra 1 is omitted when an answer is produced. ----- Steps 3 and 4 were switched because the use of the original--- algorithm for the problem 2x + y - z = 2 produces a non-minimal--- solution.  linDiophEq [2,1,-1] 2 = [[1,0,0],[0,2,0]], but the--- original algorithm produces [[1,0,0],[0,2,0],[1,1,1]].+-- When solving an inhomogeneous equation, it uses the homogeneous+-- solver after adding -b as the first element of v and by bounding+-- the first element of a to be one at each step in the computation.+-- The first element of a solution is zero if it is a solution to the+-- associated homogeneous equation, and one if it is a solution to the+-- inhomogeneous equation. -- -- The algorithm is likely to be Fortenbacher's algorithm, the one -- generalized to systems of equations by Contejean and Devie, but I--- have not been able to verified this fact.  I learned how to extend--- Contejean and Devie's results to an inhomogeneous equation by--- reading \"Effective Solutions of Linear Diophantine Equation--- Systems with an Application to Chemistry\" by David Papp and Bela--- Vizari, Rutcor Research Report RRR 28-2004, September, 2004,--- <http://rutcor.rutgers.edu/pub/rrr/reports2004/28_2004.ps>.------ The example that shows a problem with the original algorithm--- follows.  For the problem linDiophEq [2,1,-1] 2, the value of a and--- m at the beginning of the loop is:------ @---                    a                                 m---    [[0, 0, 1], [0, 1, 0], [1, 0, 0]]       []---    [[0, 1, 1], [0, 2, 0]]                  [[1, 0, 0]]---    []                                      [[1, 0, 0], [0, 2, 0]]--- @------ Consider [0, 1, 1] in a.  If you remove unnecessary branches first,--- the element will stay in a.  After performing breadth-first search,--- a will contain [1, 1, 1], which is the unwanted, non-minimal--- solution.+-- have not been able to verified this fact.  module Algebra.CommutativeMonoid.LinDiophEq (linDiophEq) where @@ -113,14 +89,32 @@  -- | The 'linDiophEq' function takes a list of integers that specifies -- the coefficients of linear Diophantine equation and a constant,--- and returns the equation's minimal, non-negative solutions.  When--- solving an inhomogeneous equation, solve the related homogeneous--- equation and add in those solutions.+-- and returns the equation's minimal, non-negative solutions.+--+-- When solving an inhomogeneous equation, the first element of a+-- solution is zero if it solves the associated homogeneous equation,+-- and one otherwise.+--+-- Thus to solve 2x + y - z = 2, use+--+-- @+-- linDiophEq [2,1,-1] 2 = [[0,0,1,1],[1,1,0,0],[0,1,0,2],[1,0,2,0]]+-- @+--+-- The two minimal solutions to the homogeneous equation are [0,1,1]+-- and [1,0,2], so any linear combinations of these solutions+-- contributes to a solution.  The solution that corresponds to+-- [1,0,0] is x = w + 1, y = v, and z = v + 2w.  The solution that+-- corresponds to [0,2,0] is x = w, y = v + 2, and z = v + 2w.+ linDiophEq :: [Int] -> Int -> [[Int]] linDiophEq [] _ = []-linDiophEq v c =-    newMinimalResults (vector n v) c (basis n) S.empty+linDiophEq v 0 =+    newMinimalResults True (vector n v) (basis n) S.empty     where n = length v+linDiophEq v c =+    newMinimalResults False (vector n (negate c:v)) (basis n) S.empty+    where n = 1 + length v  -- Construct the basis vectors for an n-dimensional space basis :: Int -> Set (Vector Int)@@ -131,35 +125,36 @@ -- This is the main loop.  -- Add elements of a that solve the equation to m and the output-newMinimalResults :: Vector Int -> Int -> Set (Vector Int) ->+-- Variable hom is true when solving a homogeneous equation+newMinimalResults :: Bool -> Vector Int -> Set (Vector Int) ->                      Set (Vector Int) -> [[Int]] newMinimalResults _ _ a _ | S.null a = []-newMinimalResults v c a m =+newMinimalResults hom v a m =     loop m (S.toList a)         -- Test each element in a     where       loop m [] =               -- When done, prepare for next iteration-          let a' = breadthFirstSearch v c a     -- Step 3-              a'' = unnecessaryBranches a' m in -- Step 4--- The original algorithm reverses these two steps.---          let a' = unnecessaryBranches a m---              a'' = breadthFirstSearch v c a' in-          newMinimalResults v c a'' m+          let a' = unnecessaryBranches a m+              a'' = breadthFirstSearch hom v a' in+          newMinimalResults hom v a'' m       loop m (x:xs)-           | prod v x == c && S.notMember x m =+           | prod v x == 0 && S.notMember x m =                elems x:loop (S.insert x m) xs -- Answer found            | otherwise =                loop m xs  -- Breadth-first search using the algorithm of Contejean and Devie-breadthFirstSearch :: Vector Int -> Int -> Set (Vector Int) -> Set (Vector Int)-breadthFirstSearch v c a =+-- Variable hom is true when solving a homogeneous equation+breadthFirstSearch :: Bool -> Vector Int ->+                      Set (Vector Int) -> Set (Vector Int)+breadthFirstSearch hom v a =     S.fold f S.empty a     where       f x acc =           foldl (flip S.insert) acc             [ x // [(k, x!k + 1)] |-              k <- indices x,-              (prod v x - c) * v!k < 0 ] -- Fortenbacher contribution+              k <- indices x,   -- When not hom, bound first element+              hom || k > 0 || x!k == 0, -- of x to be no more than one+              prod v x * v!k < 0 ] -- Fortenbacher contribution  -- Inner product prod :: Vector Int -> Vector Int -> Int