diff --git a/lagrangian.cabal b/lagrangian.cabal --- a/lagrangian.cabal +++ b/lagrangian.cabal @@ -10,25 +10,38 @@ -- PVP summary: +-+------- breaking API changes -- | | +----- non-breaking API additions -- | | | +--- code changes with no API change -version: 0.3.0.1 +version: 0.4.0.0 -- A short (one-line) description of the package. synopsis: Solve lagrange multiplier problems -- A longer description of the package. description: - Numerically solve convex lagrange multiplier problems with conjugate gradient descent. + Numerically solve convex Lagrange multiplier problems with conjugate gradient descent. . + For some background on the method of Lagrange multipliers checkout the wikipedia page + + . + Here is an example from the Wikipedia page on Lagrange multipliers + Maximize f(x, y) = x + y, subject to the constraint x^2 + y^2 = 1 + . + @ + \> maximize 0.00001 (\[x, y] -> x + y) [(\[x, y] -> x^2 + y^2) <=> 1] 2 + Right ([0.707,0.707], [-0.707]) + @ + . + For more information look here: + . For example, find the maximum entropy with the constraint that the probabilities sum to one. . @ - \> solve 0.00001 (negate . sum . map (\x -> x * log x)) [sum \<=\> 1] 3 + \> maximize 0.00001 (negate . sum . map (\\x -> x * log x)) [sum \<=\> 1] 3 Right ([0.33, 0.33, 0.33], [-0.09]) @ . The first elements of the result pair are the arguments for the - objective function at the minimum. The second elements are the lagrange multipliers. + objective function at the maximum. The second elements are the Lagrange multipliers. . -- URL for the project homepage or repository. homepage: http://github.com/jfischoff/lagrangian diff --git a/src/Numeric/AD/Lagrangian.hs b/src/Numeric/AD/Lagrangian.hs --- a/src/Numeric/AD/Lagrangian.hs +++ b/src/Numeric/AD/Lagrangian.hs @@ -3,7 +3,7 @@ -- Here is an example from the Wikipedia page on Lagrange multipliers. -- Maximize f(x, y) = x + y, subject to the constraint x^2 + y^2 = 1 -- --- >>> solve 0.00001 (\[x, y] -> x + y) [\[x, y] -> x^2 + y^2 <=> 1] 2 +-- >>> maximize 0.00001 (\[x, y] -> x + y) [(\[x, y] -> x^2 + y^2) <=> 1] 2 -- Right ([0.707,0.707], [-0.707]) -- -- The first elements of the result pair are the arguments for the objective function at the minimum. @@ -14,7 +14,8 @@ (<=>), Constraint, -- ** Solver - solve, + maximize, + minimize, -- *** Experimental features feasible) where -import Numeric.AD.Lagrangian.Internal (AD2, (<=>), solve, feasible, Constraint) +import Numeric.AD.Lagrangian.Internal (AD2, (<=>), maximize, minimize, feasible, Constraint) diff --git a/src/Numeric/AD/Lagrangian/Internal.hs b/src/Numeric/AD/Lagrangian/Internal.hs --- a/src/Numeric/AD/Lagrangian/Internal.hs +++ b/src/Numeric/AD/Lagrangian/Internal.hs @@ -27,7 +27,7 @@ -- | This is the lagrangian multiplier solver. It is assumed that the -- objective function and all of the constraints take in the -- same amount of arguments. -solve :: Double +minimize :: Double -> (forall s r. (Mode s, Mode r) => [AD2 s r Double] -> AD2 s r Double) -- ^ The function to minimize -> (forall s r. (Mode s, Mode r) => [Constraint (AD2 s r Double)] ) @@ -39,7 +39,7 @@ -> Either (Result, Statistics) (S.Vector Double, S.Vector Double) -- ^ Either an explanation of why the gradient descent failed or a pair -- containing the arguments at the minimum and the lagrange multipliers -solve tolerance toMin constraints argCount = result where +minimize tolerance toMin constraints argCount = result where -- The function to minimize for the langrangian is the squared gradient obj argsAndLams = squaredGrad (lagrangian toMin constraints argCount) argsAndLams @@ -61,6 +61,22 @@ (vs, ToleranceStatisfied, _) -> Right (S.take argCount vs, S.drop argCount vs) (_, x, y) -> Left (x, y) + +-- | Finding the maximum is the same as the minimum with the objective function inverted +maximize :: Double + -> (forall s r. (Mode s, Mode r) => [AD2 s r Double] -> AD2 s r Double) + -- ^ The function to maximize + -> (forall s r. (Mode s, Mode r) => [Constraint (AD2 s r Double)] ) + -- ^ The constraints as pairs @g \<=\> c@ which represent equations + -- of the form @g(x, y, ...) = c@ + -> Int + -- ^ The arity of the objective function which should equal the arity of + -- the constraints. + -> Either (Result, Statistics) (S.Vector Double, S.Vector Double) + -- ^ Either an explanation of why the gradient descent failed or a pair + -- containing the arguments at the minimum and the lagrange multipliers +maximize tolerance toMax constraints argCount = + minimize tolerance (negate1 . toMax) constraints argCount lagrangian :: Num a => ([a] -> a) @@ -72,7 +88,7 @@ args = take argCount argsAndLams lams = drop argCount argsAndLams - -- (g, c) <=> g(x, y, ...) = c <=> g(x, y, ...) - c = 0 + -- g(x, y, ...) = c <=> g(x, y, ...) - c = 0 appliedConstraints = fmap (\(f, c) -> f args - c) constraints -- L(x, y, ..., lam0, ...) = f(x, y, ...) + lam0 * (g0 - c0) ... diff --git a/tests/Main.hs b/tests/Main.hs --- a/tests/Main.hs +++ b/tests/Main.hs @@ -16,7 +16,7 @@ noConstraints = (fst <$> actual) @?= Right expected where - actual = solve 0.00001 f [] 1 + actual = minimize 0.00001 f [] 1 expected = S.fromList [1] f [x] = -(x - 1) ^2 @@ -24,7 +24,7 @@ -- x =~= y :: a -> a -> Bool entropyTest = (S.sum . S.map abs $ S.zipWith (-) actual expected) < 0.02 @?= True where - Right actual = fst <$> solve 0.00001 f [(\xs -> sum xs, 1)] 3 + Right actual = fst <$> maximize 0.00001 f [(\xs -> sum xs, 1)] 3 expected = S.fromList [0.33, 0.33, 0.33] f :: Floating a => [a] -> a f = negate . sum . map (\x -> x * log x)