|
10
|
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
|
