-
-
Notifications
You must be signed in to change notification settings - Fork 50
/
User_Type.jl
53 lines (41 loc) · 1.26 KB
/
User_Type.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
#=
This is an example of how to use a custom user type of with ODE module.
An additional example can be found in `test/inverface-tests.jl`
The ODE considered in this example is the simple oscilator
a' = -b
b' = a
with initial condtion a(t=0) = 0 and b(t=0) = 0.1.
=#
using ODE
using PyPlot
using Compat
typealias DataFloat Float64
type Vec2
a::DataFloat
b::DataFloat
end
# This is the minimal set of function required on the type inorder to work with
# the ODE module
@compat Base.:/(x::Vec2, y::Real) = Vec2(x.a/y, x.b/y)
@compat Base.:*(y::Real, x::Vec2) = Vec2(y * x.a, y * x.b)
@compat Base.:*(x::Vec2, y::Real) = y*x
@compat Base.:.*(y::Real, x::Vec2) = y*x
@compat Base.:+(x::Vec2, y::Vec2) = Vec2(x.a + y.a, x.b + y.b)
@compat Base.:-(x::Vec2, y::Vec2) = Vec2(x.a - y.a, x.b - y.b)
Base.norm(x::Vec2) = sqrt(x.a^2 + x.b^2)
Base.zero(x::Type{Vec2}) = Vec2(zero(DataFloat), zero(DataFloat))
ODE.isoutofdomain(x::Vec2) = isnan(x.a) || isnan(x.b)
# RHS function
f(t,y) = Vec2(-y.b, y.a)
# Initial condtions
start = Vec2(0.0, 0.1)
# Time vector going from 0 to 2*PI in 0.01 increments
time = 0:0.1:4*pi
# Solve the ODE
t, y = ode45(f, start, time)
# Plot the solution
a = map(y -> y.a, y)
b = map(y -> y.b, y)
plot(t, a, label="a(t)")
plot(t, b, label="b(t)")
legend()