Mercurial > repos > public > sbplib_julia

comparison TimeStepper.jl @ 10:bed51234616b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
add timestepper
author Ylva Rydin <ylva.rydin@telia.com>
date Mon, 17 Dec 2018 14:42:11 +0100
parents
children c4bc1165fdf7
comparison
equal deleted inserted replaced
0:b714e341a0ba 10:bed51234616b
1 abstract type TimeStepper end
2
3 # Returns v and t
4 function getState(ts::TimeStepper)
5 error("not implemented")
6 end
7
8
9 function step(ts::TimeStepper)
10 error("not implemented")
11 end
12
13 function stepN(ts::TimeStepper,N::Int)
14 for i ∈ 1:N
15 ts.step()
16 end
17 end
18
19 function stepTo(ts::TimeStepper)
20 error("Not yet implemented")
21 end
22
23 function evolve(ts::TimeStepper)
24 error("Not yet implemented")
25 end
26
27
28 mutable struct Rk4 <: TimeStepper
29 F::Function
30 k::Real
31 v::Vector
32 t::Real
33 n::UInt
34
35 function Rk4(F::Function,k::Real,v0::Vector,t0::Real)
36 # TODO: Check that F has two inputs and one output
37 v = v0
38 t = t0
39 n = 0
40 return new(F,k,v,t,n)
41 end
42 end
43
44 function getState(ts::Rk4)
45 return ts.t, ts.v
46 end
47
48 function step(ts::Rk4)
49 k1 = ts.F(ts.v,ts.t)
50 k2 = ts.F(ts.v+0.5*ts.k*k1,ts.t+0.5*ts.k)
51 k3 = ts.F(ts.v+0.5*ts.k*k2,ts.t+0.5*ts.k)
52 k4 = ts.F(ts.v+ ts.k*k3,ts.t+ ts.k)
53 ts.v = ts.v + (1/6)*(k1+2*(k2+k3)+k4)*ts.k
54
55 ts.n = ts.n + 1
56 ts.t = ts.t + ts.k
57
58 return nothing
59 end
60