@@ -7,6 +7,7 @@ using Test
7
7
@variables x (t)
8
8
D = Differential (t)
9
9
tspan = (0. , 1. )
10
+ timepoint = [0. , 0.1 , 0.2 , 0.3 , 0.4 , 0.5 , 0.6 , 0.7 , 0.8 , 0.9 , 1. ]
10
11
11
12
function expect (t, α)
12
13
return (3 / 2 * t^ (α/ 2 ) - t^ 4 )^ 2
@@ -18,7 +19,7 @@ eqs += (gamma(9)*t^(8 - α)/gamma(9 - α)) + (3/2*t^(α/2)-t^4)^3 - x^(3/2)
18
19
sys = fractional_to_ordinary (eqs, x, α, 10 ^- 7 , 1 )
19
20
20
21
prob = ODEProblem (sys, [], tspan)
21
- sol = solve (prob, radau5 (), abstol = 1e-10 , reltol = 1e-10 )
22
+ sol = solve (prob, radau5 (), saveat = timepoint, abstol = 1e-10 , reltol = 1e-10 )
22
23
23
24
time = 0
24
25
while (time <= 1 )
@@ -32,11 +33,11 @@ eqs += (gamma(9)*t^(8 - α)/gamma(9 - α)) + (3/2*t^(α/2)-t^4)^3 - x^(3/2)
32
33
sys = fractional_to_ordinary (eqs, x, α, 10 ^- 7 , 1 )
33
34
34
35
prob = ODEProblem (sys, [], tspan)
35
- sol = solve (prob, radau5 (), abstol = 1e-10 , reltol = 1e-10 )
36
+ sol = solve (prob, radau5 (), saveat = timepoint, abstol = 1e-10 , reltol = 1e-10 )
36
37
37
38
time = 0
38
39
while (time <= 1 )
39
- @test isapprox (expect (time, α), sol (time, idxs= x), atol= 1e-3 )
40
+ @test isapprox (expect (time, α), sol (time, idxs= x), atol= 1e-7 )
40
41
time += 0.1
41
42
end
42
43
@@ -46,11 +47,11 @@ eqs += (gamma(9)*t^(8 - α)/gamma(9 - α)) + (3/2*t^(α/2)-t^4)^3 - x^(3/2)
46
47
sys = fractional_to_ordinary (eqs, x, α, 10 ^- 7 , 1 )
47
48
48
49
prob = ODEProblem (sys, [], tspan)
49
- sol = solve (prob, radau5 (), abstol = 1e-10 , reltol = 1e-10 )
50
+ sol = solve (prob, radau5 (), saveat = timepoint, abstol = 1e-10 , reltol = 1e-10 )
50
51
51
52
time = 0
52
53
while (time <= 1 )
53
- @test isapprox (expect (time, α), sol (time, idxs= x), atol= 1e-3 )
54
+ @test isapprox (expect (time, α), sol (time, idxs= x), atol= 1e-7 )
54
55
time += 0.1
55
56
end
56
57
60
61
D = Differential (t)
61
62
tspan = (0. , 220. )
62
63
63
- sys = fractional_to_ordinary ([1 - 4 * x + x^ 2 * y, 3 * x - x^ 2 * y], [x, y], [1.3 , 0.8 ], 10 ^- 6 , 220 ; initials= [[1.2 , 1 ], 2.8 ])
64
+ sys = fractional_to_ordinary ([1 - 4 * x + x^ 2 * y, 3 * x - x^ 2 * y], [x, y], [1.3 , 0.8 ], 10 ^- 8 , 220 ; initials= [[1.2 , 1 ], 2.8 ])
64
65
prob = ODEProblem (sys, [], tspan)
65
66
sol = solve (prob, radau5 (), abstol = 1e-8 , reltol = 1e-8 )
66
67
67
- @test isapprox (1.0097684171 , sol (220 , idxs= x), atol= 1e-3 )
68
- @test isapprox (2.1581264031 , sol (220 , idxs= y), atol= 1e-3 )
68
+ @test isapprox (1.0097684171 , sol (220 , idxs= x), atol= 1e-5 )
69
+ @test isapprox (2.1581264031 , sol (220 , idxs= y), atol= 1e-5 )
70
+
71
+ # Testing for example 3 of Section 7
72
+ @independent_variables t
73
+ @variables x_0 (t)
74
+ D = Differential (t)
75
+ tspan = (0. , 5000. )
76
+
77
+ function expect (t)
78
+ return sqrt (2 ) * sin (t + pi / 4 )
79
+ end
80
+
81
+ sys = linear_fractional_to_ordinary ([3 , 2.5 , 2 , 1 , .5 , 0 ], [1 , 1 , 1 , 4 , 1 , 4 ], 6 * cos (t), 10 ^- 5 , 5000 ; initials= [1 , 1 , - 1 ])
82
+ prob = ODEProblem (sys, [], tspan)
83
+ sol = solve (prob, radau5 (), abstol = 1e-5 , reltol = 1e-5 )
84
+
85
+ @test isapprox (expect (5000 ), sol (5000 , idxs= x_0), atol= 1e-6 )
0 commit comments