jfischoff / lagrangian

Solve lagrangian multiplier problems
BSD 3-Clause "New" or "Revised" License
5 stars 4 forks source link

Increase type safety #4

Open jfischoff opened 9 years ago

jfischoff commented 9 years ago

I really find the interface annoying.

I initially blew off the need to explicitly specify the arity of the function in minimize but caused a demo to crash when presenting the maxent library, so, I think.

I think it is possible to come up with a safer interface.

jfischoff commented 9 years ago

I'm thinking the following interface:

minimize :: forall (argCount :: Nat) (constraintCount :: Nat). 
            ( ( forall a.
              , Floating a
              , Traversable (f n)
              ) 
            => f argCount a -> a
            )
         -- ^ The objective function to minimize
         -> Vec constraintCount Constraint
         -- ^ A list of constraints @g \<=\> c@ corresponding to equations of
         -- the form @g(x, y, ...) = c@
         -> Double
         -- ^ Stop iterating when the largest component of the gradient is
         -- smaller than this value
         -> Either 
              (Result, Statistics) 
              (V.Vector Double, V.Vector Double)

Anyway the idea is to sized indexed types to get the type safety.

I'm not sure I will be able to generalize it for the traversable, hopefully it type checks.

jfischoff commented 9 years ago

Actually I think there should be a safe wrapper around the current interface.

ericpashman commented 9 years ago

I think a wrapper makes sense. A type-safe interface is going to put an extra burden on the user, since it won't be possible to pass in bog-standard objective functions like \[x, y, z] -> x^2 + y - z. Meanwhile we should probably try to make the existing approach a bit safer, e.g., by doing what we can do verify that the objective function can actually handle the specified number of arguments.

jfischoff commented 9 years ago

Added a first pass in d6667d8ad500e1766dd401a1c06fbd77a69e3514 haven't tested it

jfischoff commented 9 years ago

Whoops forgot to checkin the new file 2928ec904e4179aeae0c9a2f26533b58c31982bd

ericpashman commented 9 years ago

Good work.

But the forall f ... . Traversable (f argCount), ... => ... part of the signature on the objective function is a problem. How can you write a function with that type? I don't think it works to move the f and the Traversable constraint out from behind the rank-2 quantifier; f will have to unify with IList regardless since the (value-level) function definition uses the IList constructor.

So the easy fix is to use forall a. Floating a => IList a -> a instead. Then you can write objective functions like this

>> let g (IList xs) = negate . sum . map (\x -> x * log x) $ xs

and call the optimizers like this

>> import Numeric.AD.Lagrangian.Safe as S
>> :set -XDataKinds
>> S.minimize (g :: Floating a => IList 5 a -> a) [sum <=> 1] (1e-10)
Right (fromList [0.19999999999818122,0.19999999999818122,0.19999999999818122,0.19999999999818122,0.19999999999818122],fromList [-0.6094379124449681])

That works. Is it more or less what you had in mind?

I'll submit a patch as a pull request so you can see whether there's a more general way to do things before merging.

ericpashman commented 9 years ago

Oh, and IList should probably be a newtype. Having made it one, you can also use GHC's snazzy deriving features to write the necessary instances.

Edit: I updated the pull request above to make this change.

jfischoff commented 9 years ago

I think your are right about the higher rank type.

I was trying to not force people to use a particular sized indexed data structure since there are many possible ones. I'm going to just have this work with Vec n a.

The IList thing is silly.

ericpashman commented 9 years ago

Yeah, I guess what we most want is to prevent users from passing in objective functions that only work for inputs of a fixed dimension and then mis-specifying the dimension and watching everything blow up at run-time, right?

There are ways other than type-level programming to do that, particularly since we already have that nasty unsafePerformIO in there. For instance, we could wrap the call to optimize in some exception-handling logic that performs the optimization only if the objective function works with the specified arity, returning a Left with some sort of error otherwise:

>> import Control.Exception
>> let f [x] = x^2
>>  try (return $! f [1,2]) :: IO (Either SomeException Int)
Left <interactive>:2:5-13: Non-exhaustive patterns in function f
>> -- Broken, return a `Left` with an error.
>> try (return $! f [1]) :: IO (Either SomeException Int)
Right 1
>> -- Works, perform the optimization.

That sort of solution seems god-awful, but it would lend some safety even to the unsafe interface.

Edit: That wasn't very clear. I mean that we could do something like try (return $! f (replicate argCount 0)) :: IO (Either SomeException Int) as a test. So if argCount gives an arity the function can deal with, the try returns a Right and you go ahead and perform the optimization; if not, it returns a Left and you abort.

jfischoff commented 9 years ago

Yeah, I'm realizing this needs some more thought. I'm going to put this into another branch.

I want to get out the changes that you have made first.

Then I want to redo the linear part of maxent to not use [[a]] but an unboxed Vector or Repa or Accelerate: something efficient.