From 59c49750fd2455574fe4e67ddd7e67a20321c8a8 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Wed, 20 Feb 2013 18:56:13 -0500 Subject: [PATCH] Use RealField/RealRing where possible instead of their constituents. Convert src/ODE/IVP.hs to numeric-prelude. Add RealRing/RealField instanced to BigFloat. --- src/BigFloat.hs | 7 +++++++ src/ODE/IVP.hs | 20 +++++++++++++++----- src/Roots/Fast.hs | 29 +++++++++++++---------------- src/Roots/Simple.hs | 43 ++++++++++++++----------------------------- 4 files changed, 49 insertions(+), 50 deletions(-) diff --git a/src/BigFloat.hs b/src/BigFloat.hs index e23ebe6..0bebab9 100644 --- a/src/BigFloat.hs +++ b/src/BigFloat.hs @@ -10,6 +10,8 @@ import NumericPrelude hiding (abs) import qualified Algebra.Absolute as Absolute import qualified Algebra.Additive as Additive import qualified Algebra.Field as Field +import qualified Algebra.RealField as RealField +import qualified Algebra.RealRing as RealRing import qualified Algebra.Ring as Ring import qualified Algebra.ToRational as ToRational import qualified Algebra.ZeroTestable as ZeroTestable @@ -38,3 +40,8 @@ instance Epsilon e => ZeroTestable.C (BigFloat e) where instance Epsilon e => ToRational.C (BigFloat e) where toRational = fromRational . P.toRational + +instance Epsilon e => RealRing.C (BigFloat e) where + floor = fromInteger . P.floor + +instance Epsilon e => RealField.C (BigFloat e) diff --git a/src/ODE/IVP.hs b/src/ODE/IVP.hs index 6bde763..b4fd3dc 100644 --- a/src/ODE/IVP.hs +++ b/src/ODE/IVP.hs @@ -1,4 +1,5 @@ {-# LANGUAGE ScopedTypeVariables #-} +{-# LANGUAGE RebindableSyntax #-} -- | Numerical solution of the initial value problem, -- @@ -12,7 +13,12 @@ module ODE.IVP where import Misc (partition) - +import NumericPrelude hiding (abs) +import Algebra.Absolute (abs) +import qualified Algebra.Field as Field +import qualified Algebra.ToInteger as ToInteger +import qualified Algebra.ToRational as ToRational +import qualified Algebra.RealField as RealField -- | A single iteration of Euler's method over the interval -- [$x0$, $x0$+$h$]. @@ -26,7 +32,7 @@ import Misc (partition) -- >>> eulers_method1 x0 y0 f h -- 2.0 -- -eulers_method1 :: (RealFrac a, RealFrac b) +eulers_method1 :: (Field.C a, ToRational.C a, Field.C b) => a -- ^ x0, the initial point -> b -- ^ y0, the initial value at x0 -> (a -> b -> b) -- ^ The function f(x,y) @@ -35,7 +41,7 @@ eulers_method1 :: (RealFrac a, RealFrac b) eulers_method1 x0 y0 f h = y0 + h'*y' where - h' = realToFrac h + h' = fromRational'$ toRational h y' = (f x0 y0) @@ -59,7 +65,11 @@ eulers_method1 x0 y0 f h = -- >>> head ys == y0 -- True -- -eulers_method :: forall a b c. (RealFrac a, RealFrac b, Integral c) +eulers_method :: forall a b c. (Field.C a, + ToRational.C a, + Field.C b, + ToInteger.C c, + Enum c) => a -- ^ x0, the initial point -> b -- ^ y0, the initial value at x0 -> a -- ^ xN, the terminal point @@ -102,7 +112,7 @@ eulers_method x0 y0 xN f n = -- >>> ys == ys' -- True -- -eulers_methodH :: (RealFrac a, RealFrac b) +eulers_methodH :: (RealField.C a, ToRational.C a, Field.C b) => a -- ^ x0, the initial point -> b -- ^ y0, the initial value at x0 -> a -- ^ xN, the terminal point diff --git a/src/Roots/Fast.hs b/src/Roots/Fast.hs index 47fa512..78c299a 100644 --- a/src/Roots/Fast.hs +++ b/src/Roots/Fast.hs @@ -13,15 +13,14 @@ import Data.List (find) import Normed import NumericPrelude hiding (abs) -import Algebra.Absolute -import Algebra.Field -import Algebra.Ring - -has_root :: (Algebra.Field.C a, - Ord a, - Algebra.Ring.C b, - Ord b, - Algebra.Absolute.C b) +import qualified Algebra.Absolute as Absolute +import qualified Algebra.Field as Field +import qualified Algebra.RealRing as RealRing +import qualified Algebra.RealField as RealField + +has_root :: (RealField.C a, + RealRing.C b, + Absolute.C b) => (a -> b) -- ^ The function @f@ -> a -- ^ The \"left\" endpoint, @a@ -> a -- ^ The \"right\" endpoint, @b@ @@ -61,11 +60,9 @@ has_root f a b epsilon f_of_a f_of_b = c = (a + b)/2 -bisect :: (Algebra.Field.C a, - Ord a, - Algebra.Ring.C b, - Ord b, - Algebra.Absolute.C b) +bisect :: (RealField.C a, + RealRing.C b, + Absolute.C b) => (a -> b) -- ^ The function @f@ whose root we seek -> a -- ^ The \"left\" endpoint of the interval, @a@ -> a -- ^ The \"right\" endpoint of the interval, @b@ @@ -119,8 +116,8 @@ fixed_point_iterations f x0 = -- We also return the number of iterations required. -- fixed_point_with_iterations :: (Normed a, - Algebra.Field.C b, - Algebra.Absolute.C b, + Field.C b, + Absolute.C b, Ord b) => (a -> a) -- ^ The function @f@ to iterate. -> b -- ^ The tolerance, @epsilon@. diff --git a/src/Roots/Simple.hs b/src/Roots/Simple.hs index 44d3d62..0a1debf 100644 --- a/src/Roots/Simple.hs +++ b/src/Roots/Simple.hs @@ -18,9 +18,11 @@ import Normed import qualified Roots.Fast as F import NumericPrelude hiding (abs) -import Algebra.Absolute -import Algebra.Field -import Algebra.Ring +import qualified Algebra.Absolute as Absolute +import Algebra.Absolute (abs) +import qualified Algebra.Field as Field +import qualified Algebra.RealField as RealField +import qualified Algebra.RealRing as RealRing -- | Does the (continuous) function @f@ have a root on the interval -- [a,b]? If f(a) <] 0 and f(b) ]> 0, we know that there's a root in @@ -41,11 +43,7 @@ import Algebra.Ring -- >>> has_root cos (-2) 2 (Just 0.001) -- True -- -has_root :: (Algebra.Field.C a, - Ord a, - Algebra.Ring.C b, - Algebra.Absolute.C b, - Ord b) +has_root :: (RealField.C a, RealRing.C b) => (a -> b) -- ^ The function @f@ -> a -- ^ The \"left\" endpoint, @a@ -> a -- ^ The \"right\" endpoint, @b@ @@ -75,11 +73,7 @@ has_root f a b epsilon = -- >>> bisect sin (-1) 1 0.001 -- Just 0.0 -- -bisect :: (Algebra.Field.C a, - Ord a, - Algebra.Ring.C b, - Algebra.Absolute.C b, - Ord b) +bisect :: (RealField.C a, RealRing.C b) => (a -> b) -- ^ The function @f@ whose root we seek -> a -- ^ The \"left\" endpoint of the interval, @a@ -> a -- ^ The \"right\" endpoint of the interval, @b@ @@ -93,10 +87,7 @@ bisect f a b epsilon = -- at x0. We delegate to the version that returns the number of -- iterations and simply discard the number of iterations. -- -fixed_point :: (Normed a, - Algebra.Field.C b, - Algebra.Absolute.C b, - Ord b) +fixed_point :: (Normed a, RealField.C b) => (a -> a) -- ^ The function @f@ to iterate. -> b -- ^ The tolerance, @epsilon@. -> a -- ^ The initial value @x0@. @@ -109,10 +100,7 @@ fixed_point f epsilon x0 = -- the function @f@ with the search starting at x0 and tolerance -- @epsilon@. We delegate to the version that returns the number of -- iterations and simply discard the fixed point. -fixed_point_iteration_count :: (Normed a, - Algebra.Field.C b, - Algebra.Absolute.C b, - Ord b) +fixed_point_iteration_count :: (Normed a, RealField.C b) => (a -> a) -- ^ The function @f@ to iterate. -> b -- ^ The tolerance, @epsilon@. -> a -- ^ The initial value @x0@. @@ -130,10 +118,7 @@ fixed_point_iteration_count f epsilon x0 = -- -- This is used to determine the rate of convergence. -- -fixed_point_error_ratios :: (Normed a, - Algebra.Field.C b, - Algebra.Absolute.C b, - Ord b) +fixed_point_error_ratios :: (Normed a, RealField.C b) => (a -> a) -- ^ The function @f@ to iterate. -> a -- ^ The initial value @x0@. -> a -- ^ The true solution, @x_star@. @@ -160,7 +145,7 @@ fixed_point_error_ratios f x0 x_star p = -- >>> tail $ take 4 $ newton_iterations f f' 2 -- [1.6806282722513088,1.4307389882390624,1.2549709561094362] -- -newton_iterations :: (Algebra.Field.C a) +newton_iterations :: (Field.C a) => (a -> a) -- ^ The function @f@ whose root we seek -> (a -> a) -- ^ The derivative of @f@ -> a -- ^ Initial guess, x-naught @@ -195,7 +180,7 @@ newton_iterations f f' x0 = -- >>> abs (f root) < eps -- True -- -newtons_method :: (Algebra.Field.C a, Algebra.Absolute.C a, Ord a) +newtons_method :: (RealField.C a) => (a -> a) -- ^ The function @f@ whose root we seek -> (a -> a) -- ^ The derivative of @f@ -> a -- ^ The tolerance epsilon @@ -245,7 +230,7 @@ iterate2 f x0 x1 = -- >>> take 4 $ secant_iterations f 2 1 -- [2.0,1.0,1.0161290322580645,1.190577768676638] -- -secant_iterations :: (Algebra.Field.C a) +secant_iterations :: (Field.C a) => (a -> a) -- ^ The function @f@ whose root we seek -> a -- ^ Initial guess, x-naught -> a -- ^ Second initial guess, x-one @@ -273,7 +258,7 @@ secant_iterations f x0 x1 = -- >>> abs (f root) < (1/10^9) -- True -- -secant_method :: (Algebra.Field.C a, Algebra.Absolute.C a, Ord a) +secant_method :: (RealField.C a) => (a -> a) -- ^ The function @f@ whose root we seek -> a -- ^ The tolerance epsilon -> a -- ^ Initial guess, x-naught -- 2.43.2