diff options
Diffstat (limited to 'packages/hmatrix/examples/ode.hs')
-rw-r--r-- | packages/hmatrix/examples/ode.hs | 50 |
1 files changed, 0 insertions, 50 deletions
diff --git a/packages/hmatrix/examples/ode.hs b/packages/hmatrix/examples/ode.hs deleted file mode 100644 index dc6e0ec..0000000 --- a/packages/hmatrix/examples/ode.hs +++ /dev/null | |||
@@ -1,50 +0,0 @@ | |||
1 | {-# LANGUAGE ViewPatterns #-} | ||
2 | import Numeric.GSL.ODE | ||
3 | import Numeric.LinearAlgebra | ||
4 | import Graphics.Plot | ||
5 | import Debug.Trace(trace) | ||
6 | debug x = trace (show x) x | ||
7 | |||
8 | vanderpol mu = do | ||
9 | let xdot mu t [x,v] = [v, -x + mu * v * (1-x^2)] | ||
10 | ts = linspace 1000 (0,50) | ||
11 | sol = toColumns $ odeSolve (xdot mu) [1,0] ts | ||
12 | mplot (ts : sol) | ||
13 | mplot sol | ||
14 | |||
15 | |||
16 | harmonic w d = do | ||
17 | let xdot w d t [x,v] = [v, a*x + b*v] where a = -w^2; b = -2*d*w | ||
18 | ts = linspace 100 (0,20) | ||
19 | sol = odeSolve (xdot w d) [1,0] ts | ||
20 | mplot (ts : toColumns sol) | ||
21 | |||
22 | |||
23 | kepler v a = mplot (take 2 $ toColumns sol) where | ||
24 | xdot t [x,y,vx,vy] = [vx,vy,x*k,y*k] | ||
25 | where g=1 | ||
26 | k=(-g)*(x*x+y*y)**(-1.5) | ||
27 | ts = linspace 100 (0,30) | ||
28 | sol = odeSolve xdot [4, 0, v * cos (a*degree), v * sin (a*degree)] ts | ||
29 | degree = pi/180 | ||
30 | |||
31 | |||
32 | main = do | ||
33 | vanderpol 2 | ||
34 | harmonic 1 0 | ||
35 | harmonic 1 0.1 | ||
36 | kepler 0.3 60 | ||
37 | kepler 0.4 70 | ||
38 | vanderpol' 2 | ||
39 | |||
40 | -- example of odeSolveV with jacobian | ||
41 | vanderpol' mu = do | ||
42 | let xdot mu t (toList->[x,v]) = fromList [v, -x + mu * v * (1-x^2)] | ||
43 | jac t (toList->[x,v]) = (2><2) [ 0 , 1 | ||
44 | , -1-2*x*v*mu, mu*(1-x**2) ] | ||
45 | ts = linspace 1000 (0,50) | ||
46 | hi = (ts@>1 - ts@>0)/100 | ||
47 | sol = toColumns $ odeSolveV (MSBDF jac) hi 1E-8 1E-8 (xdot mu) (fromList [1,0]) ts | ||
48 | mplot sol | ||
49 | |||
50 | |||