Implement forward substitute in terms of a fold.

[numerical-analysis.git] / src / Linear / Iteration.hs

diff --git a/src/Linear/Iteration.hs b/src/Linear/Iteration.hs
index 6fb279c78bbb64fc5122070c82bce6f956e16a14..6117770a52313c21546ed57e729854a7d1fc5f8c 100644 (file)
--- a/src/Linear/Iteration.hs
+++ b/src/Linear/Iteration.hs
@@ -27,6 +27,7 @@ import qualified Algebra.RealField as RealField ( C )
 import qualified Algebra.ToRational as ToRational ( C )
 
 import Linear.Matrix (
+  Col,
   Mat(..),
   (!!!),
   (*),
@@ -41,12 +42,12 @@ import Normed ( Normed(..) )
 -- | A generalized implementation for Jacobi, Gauss-Seidel, etc. All
 --   that we really need to know is how to construct the matrix M, so we
 --   take a function that does it as an argument.
-classical_iteration :: (Eq a, Field.C a, Arity m)
-                 => (Mat m m  a -> Mat m m  a)
-                 -> Mat m m  a
-                 -> Mat m N1 a
-                 -> Mat m N1 a
-                 -> Mat m N1 a
+classical_iteration :: (Eq a, Field.C a, m ~ S n, Arity n)
+                 => (Mat m m a -> Mat m m a)
+                 -> Mat m m a
+                 -> Col m a
+                 -> Col m a
+                 -> Col m a
 classical_iteration m_function matrix b x_current =
   x_next
   where
@@ -59,8 +60,8 @@ classical_iteration m_function matrix b x_current =
 
 -- | Perform one iteration of successive over-relaxation.
 --
-sor_iteration :: forall m a.
-                 (Eq a, Field.C a, Arity m)
+sor_iteration :: forall m n a.
+                 (Eq a, Field.C a, m ~ S n, Arity n)
               => a -- ^ Omega
               -> Mat m m  a -- ^ Matrix A
               -> Mat m N1 a -- ^ Vector b
@@ -79,7 +80,7 @@ sor_iteration omega =
 
 -- | Compute an infinite list of SOR iterations starting with the
 --   vector x0.
-sor_iterations :: (Eq a, Field.C a, Arity m)
+sor_iterations :: (Eq a, Field.C a, m ~ S n, Arity n)
                => a
                -> Mat m m  a
                -> Mat m N1 a
@@ -90,7 +91,7 @@ sor_iterations omega matrix b =
 
 
 -- | Perform one iteration of Gauss-Seidel.
-gauss_seidel_iteration :: (Eq a, Field.C a, Arity m)
+gauss_seidel_iteration :: (Eq a, Field.C a, m ~ S n, Arity n)
                        => Mat m m  a
                        -> Mat m N1 a
                        -> Mat m N1 a
@@ -100,7 +101,7 @@ gauss_seidel_iteration = sor_iteration one
 
 -- | Compute an infinite list of Gauss-Seidel iterations starting with
 --   the vector x0.
-gauss_seidel_iterations :: (Eq a, Field.C a, Arity m)
+gauss_seidel_iterations :: (Eq a, Field.C a, m ~ S n, Arity n)
                         => Mat m m  a
                         -> Mat m N1 a
                         -> Mat m N1 a
@@ -126,7 +127,7 @@ gauss_seidel_iterations matrix b =
 --   >>> jacobi_iteration m b x1
 --   ((0.0),(0.25))
 --
-jacobi_iteration :: (Eq a, Field.C a, Arity m)
+jacobi_iteration :: (Eq a, Field.C a, m ~ S n, Arity n)
                  => Mat m m  a
                  -> Mat m N1 a
                  -> Mat m N1 a
@@ -137,7 +138,7 @@ jacobi_iteration =
 
 -- | Compute an infinite list of Jacobi iterations starting with the
 --   vector x0.
-jacobi_iterations :: (Eq a, Field.C a, Arity m)
+jacobi_iterations :: (Eq a, Field.C a, m ~ S n, Arity n)
                   => Mat m m  a
                   -> Mat m N1 a
                   -> Mat m N1 a