diff --git a/optimization.cabal b/optimization.cabal
--- a/optimization.cabal
+++ b/optimization.cabal
@@ -1,6 +1,6 @@
 name:          optimization
 category:      Math
-version:       0.1
+version:       0.1.1
 license:       BSD3
 cabal-version: >= 1.10
 license-file:  LICENSE
diff --git a/src/Optimization/Constrained/Penalty.hs b/src/Optimization/Constrained/Penalty.hs
--- a/src/Optimization/Constrained/Penalty.hs
+++ b/src/Optimization/Constrained/Penalty.hs
@@ -61,6 +61,7 @@
   where go x0 l0 = let l1 = fmap (*alpha) l0
                        x1 = head $ drop 100 $ minX (FU $ \x -> augLagrangian opt x (fmap auto l1)) x0
                    in x1 : go x1 l1
+{-# INLINEABLE minimize #-}
 
 -- | Maximize the given constrained optimization problem
 maximize :: (Functor f, Num a, Ord a, g ~ V)
@@ -72,12 +73,14 @@
          -> [f a]                      -- ^ Optimizing iterates
 maximize minX (Opt (FU f) gs hs) alpha =
     minimize minX (Opt (FU $ negate . f) gs hs) alpha
+{-# INLINEABLE maximize #-}
 
 -- | The Lagrangian for the given constrained optimization
 lagrangian :: (Num a) => Opt f a
            -> (forall s. Mode s => f (AD s a) -> V (AD s a) -> AD s a)
 lagrangian (Opt (FU f) gs hs) x l =
     f x - V.sum (V.zipWith (\lamb (FU g)->lamb * g x) l gs)
+{-# INLINEABLE lagrangian #-}
 
 -- | The augmented Lagrangian for the given constrained optimization
 augLagrangian :: (Num a, Ord a) => Opt f a
@@ -86,3 +89,4 @@
     f x + V.sum (V.zipWith (*) l $ V.concat [gs', hs'])
   where gs' = V.map (\(FU g) -> (g x)^2) gs
         hs' = V.map (\(FU h) -> (max 0 $ h x)^2) hs
+{-# INLINEABLE augLagrangian #-}
diff --git a/src/Optimization/LineSearch.hs b/src/Optimization/LineSearch.hs
--- a/src/Optimization/LineSearch.hs
+++ b/src/Optimization/LineSearch.hs
@@ -30,8 +30,8 @@
 import Prelude hiding (pred)
 import Linear
 
--- | A 'LineSearch' method 'search df p x' determines a step size
--- in direction 'p' from point 'x' for function 'f' with gradient 'df'
+-- | A line search method @search df p x@ determines a step size
+-- in direction @p@ from point @x@ for function @f@ with gradient @df@
 type LineSearch f a = (f a -> f a) -> f a -> f a -> a
 
 -- | Armijo condition
@@ -63,6 +63,7 @@
   where nonzero (x:xs) | not $ x > 0 = error "Backtracking search failed: alpha=0" -- FIXME
                        | otherwise   = x : nonzero xs
         nonzero [] = error "Backtracking search failed: no more iterates"
+{-# INLINEABLE backtrackingSearch #-}
 
 -- | Armijo backtracking line search algorithm
 --
@@ -73,6 +74,7 @@
              => a -> a -> a -> (f a -> a) -> LineSearch f a
 armijoSearch gamma alpha c1 f df p x =
     backtrackingSearch gamma alpha (armijo c1 f df x p) df p x
+{-# INLINEABLE armijoSearch #-}
 
 -- | Wolfe backtracking line search algorithm
 --
@@ -84,6 +86,7 @@
 wolfeSearch gamma alpha c1 c2 f df p x =
     backtrackingSearch gamma alpha wolfe df p x
   where wolfe a = armijo c1 f df p x a && curvature c2 df x p a
+{-# INLINEABLE wolfeSearch #-}
 
 -- | Line search by Newton's method
 newtonSearch :: (Num a) => LineSearch f a
@@ -98,3 +101,4 @@
 -- @constantSearch c@ always chooses a step-size @c@.
 constantSearch :: a -> LineSearch f a
 constantSearch c _ _ _ = c
+{-# INLINEABLE constantSearch #-}
diff --git a/src/Optimization/LineSearch/BFGS.hs b/src/Optimization/LineSearch/BFGS.hs
--- a/src/Optimization/LineSearch/BFGS.hs
+++ b/src/Optimization/LineSearch/BFGS.hs
@@ -29,3 +29,5 @@
                          v = i !-! rho *!! outer s y
                          b1 = u !*! b0 !*! v !+! rho *!! outer s s
                      in x1 : go b1 x1
+{-# INLINABLE bfgs #-}    
+
diff --git a/src/Optimization/LineSearch/BarzilaiBorwein.hs b/src/Optimization/LineSearch/BarzilaiBorwein.hs
--- a/src/Optimization/LineSearch/BarzilaiBorwein.hs
+++ b/src/Optimization/LineSearch/BarzilaiBorwein.hs
@@ -15,3 +15,4 @@
                    in if nearZero (z `dot` z)
                         then [x2]
                         else x2 : go x1 x2
+{-# INLINABLE barzilaiBorwein #-}  
diff --git a/src/Optimization/LineSearch/ConjugateGradient.hs b/src/Optimization/LineSearch/ConjugateGradient.hs
--- a/src/Optimization/LineSearch/ConjugateGradient.hs
+++ b/src/Optimization/LineSearch/ConjugateGradient.hs
@@ -13,9 +13,9 @@
 import Optimization.LineSearch
 import Linear
 
--- | A beta expression 'beta df0 df1 p' is an expression for the
--- conjugate direction contribution given the derivative 'df0' and
--- direction 'p' for iteration 'k', 'df1' for iteration 'k+1'
+-- | A beta expression @beta df0 df1 p@ is an expression for the
+-- conjugate direction contribution given the derivative @df0@ and
+-- direction @p@ for iteration @k@, @df1@ for iteration @k+1@
 type Beta f a = f a -> f a -> f a -> a
 
 -- | Conjugate gradient method with given beta and line search method
diff --git a/src/Optimization/LineSearch/MirrorDescent.hs b/src/Optimization/LineSearch/MirrorDescent.hs
--- a/src/Optimization/LineSearch/MirrorDescent.hs
+++ b/src/Optimization/LineSearch/MirrorDescent.hs
@@ -21,3 +21,5 @@
                     y1 = dPsi x0 ^-^ t0 *^ df x0
                     x1 = dPsiStar y1
                 in x1 : go y1
+{-# INLINEABLE mirrorDescent #-}
+
diff --git a/src/Optimization/TrustRegion/Nesterov1983.hs b/src/Optimization/TrustRegion/Nesterov1983.hs
--- a/src/Optimization/TrustRegion/Nesterov1983.hs
+++ b/src/Optimization/TrustRegion/Nesterov1983.hs
@@ -10,7 +10,6 @@
 -- gradient method, first described in 1983. This method requires
 -- knowledge of the Lipschitz constant @l@ of the gradient, the condition
 -- number @kappa@, as well as an initial step size @alpha0@ in @(0,1)@.
-{-# INLINEABLE optimalGradient #-}
 optimalGradient :: (Additive f, Functor f, Ord a, Floating a, Epsilon a)
                 => a -> a -> (f a -> f a) -> a -> f a -> [f a]
 optimalGradient kappa l df a0' x0' = go a0' x0' x0'
@@ -22,9 +21,10 @@
                           b1 = a0 * (1 - a0) / (a0^2 + a1)
                           y1 = x1 ^+^ b1 *^ (x1 ^-^ x0)
                       in x1 : go a0 x1 y1
+{-# INLINEABLE optimalGradient #-}
 
--- | 'quadratic a b c' is the real solutions to a quadratic equation
--- 'a x^2 + b x + c == 0'
+-- | @quadratic a b c@ is the real solutions to a quadratic equation
+-- @a x^2 + b x + c == 0@
 quadratic :: (Ord a, Floating a, Epsilon a)
           => a -> a -> a -> [a]
 quadratic a b c
@@ -33,3 +33,4 @@
     | otherwise      = [ (-b + sqrt discr) / 2 / a
                        , (-b - sqrt discr) / 2 / a ]
   where discr = b^2 - 4*a*c
+{-# INLINEABLE quadratic #-}
diff --git a/src/Optimization/TrustRegion/Newton.hs b/src/Optimization/TrustRegion/Newton.hs
--- a/src/Optimization/TrustRegion/Newton.hs
+++ b/src/Optimization/TrustRegion/Newton.hs
@@ -13,11 +13,11 @@
 import Linear
 
 -- | Newton's method
-{-# INLINEABLE newton #-}
 newton :: (Num a, Ord a, Additive f, Metric f, Foldable f)
        => (f a -> f a) -> (f a -> f (f a)) -> f a -> [f a]
 newton df ddfInv x0 = iterate go x0
   where go x = x ^-^ ddfInv x !* df x
+{-# INLINEABLE newton #-}
 
 -- | Inverse by iterative method of Ben-Israel and Cohen
 -- with given starting condition
@@ -26,12 +26,15 @@
         => m (m a) -> m (m a) -> [m (m a)]
 bicInv' a0 a = iterate go a0
   where go ak = 2 *!! ak !-! ak !*! a !*! ak
+{-# INLINEABLE bicInv' #-}
 
 -- | Inverse by iterative method of Ben-Israel and Cohen
--- starting from 'alpha A^T'. Alpha should be set such that
+-- starting from @alpha A^T@. Alpha should be set such that
 -- 0 < alpha < 2/sigma^2 where sigma denotes the largest singular
 -- value of A
 bicInv :: (Functor m, Distributive m, Additive m,
            Applicative m, Apply m, Foldable m, Conjugate a)
        => a -> m (m a) -> [m (m a)]
 bicInv alpha a = bicInv' (alpha *!! adjoint a) a
+{-# INLINEABLE bicInv #-}
+