diff options
Diffstat (limited to 'packages/gsl/src/Numeric/GSL/Polynomials.hs')
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Polynomials.hs | 57 |
1 files changed, 57 insertions, 0 deletions
diff --git a/packages/gsl/src/Numeric/GSL/Polynomials.hs b/packages/gsl/src/Numeric/GSL/Polynomials.hs new file mode 100644 index 0000000..b1be85d --- /dev/null +++ b/packages/gsl/src/Numeric/GSL/Polynomials.hs | |||
@@ -0,0 +1,57 @@ | |||
1 | {- | | ||
2 | Module : Numeric.GSL.Polynomials | ||
3 | Copyright : (c) Alberto Ruiz 2006 | ||
4 | License : GPL | ||
5 | Maintainer : Alberto Ruiz | ||
6 | Stability : provisional | ||
7 | |||
8 | Polynomials. | ||
9 | |||
10 | <http://www.gnu.org/software/gsl/manual/html_node/General-Polynomial-Equations.html#General-Polynomial-Equations> | ||
11 | |||
12 | -} | ||
13 | |||
14 | |||
15 | module Numeric.GSL.Polynomials ( | ||
16 | polySolve | ||
17 | ) where | ||
18 | |||
19 | import Data.Packed | ||
20 | import Numeric.GSL.Internal | ||
21 | import Data.Complex | ||
22 | import System.IO.Unsafe (unsafePerformIO) | ||
23 | |||
24 | #if __GLASGOW_HASKELL__ >= 704 | ||
25 | import Foreign.C.Types (CInt(..)) | ||
26 | #endif | ||
27 | |||
28 | {- | Solution of general polynomial equations, using /gsl_poly_complex_solve/. | ||
29 | |||
30 | For example, the three solutions of x^3 + 8 = 0 | ||
31 | |||
32 | >>> polySolve [8,0,0,1] | ||
33 | [(-2.0) :+ 0.0,1.0 :+ 1.7320508075688776,1.0 :+ (-1.7320508075688776)] | ||
34 | |||
35 | |||
36 | The example in the GSL manual: To find the roots of x^5 -1 = 0: | ||
37 | |||
38 | >>> polySolve [-1, 0, 0, 0, 0, 1] | ||
39 | [(-0.8090169943749472) :+ 0.5877852522924731, | ||
40 | (-0.8090169943749472) :+ (-0.5877852522924731), | ||
41 | 0.30901699437494756 :+ 0.9510565162951535, | ||
42 | 0.30901699437494756 :+ (-0.9510565162951535), | ||
43 | 1.0000000000000002 :+ 0.0] | ||
44 | |||
45 | -} | ||
46 | polySolve :: [Double] -> [Complex Double] | ||
47 | polySolve = toList . polySolve' . fromList | ||
48 | |||
49 | polySolve' :: Vector Double -> Vector (Complex Double) | ||
50 | polySolve' v | dim v > 1 = unsafePerformIO $ do | ||
51 | r <- createVector (dim v-1) | ||
52 | app2 c_polySolve vec v vec r "polySolve" | ||
53 | return r | ||
54 | | otherwise = error "polySolve on a polynomial of degree zero" | ||
55 | |||
56 | foreign import ccall unsafe "gsl-aux.h polySolve" c_polySolve:: TV (TCV Res) | ||
57 | |||