Skip to content

Commit

Permalink
Fixing mistakes related to the constants a61 to a64
Browse files Browse the repository at this point in the history
Fixing error in perform step. This has greatly increased the order.

Update ode_convergence_tests.jl
Remove DS_Store
  • Loading branch information
ChrisRackauckas committed Jul 11, 2018
1 parent 8606aa3 commit 12fcc1b
Show file tree
Hide file tree
Showing 3 changed files with 20 additions and 14 deletions.
14 changes: 10 additions & 4 deletions src/perform_step/fixed_timestep_perform_step.jl
Original file line number Diff line number Diff line change
Expand Up @@ -414,11 +414,14 @@ end
w = integrator.alg.w
v = w*dt
## Formula by Z.A. Anastassi, see the Anas5 caches in tableaus/low_order_rk_tableaus.jl for the full citation.
a65 = (-8000//1071)*(-a43*v^5 + 6*tan(v)*v^4 + 24*v^3 - 72*tan(v)*v^2 - 144*v + 144*tan(v))/(v^5*(a43*tan(v)*v + 12 - 10*a43))
a65 = (-8000//1071)*(-a43*(v^5) + 6*tan(v)*(v^4) + 24*(v^3) - 72*tan(v)*(v^2) - 144*v + 144*tan(v))/((v^5)*(a43*tan(v)*v + 12 - 10*a43))
a61 += (-119//200)*a65
a63 += (189//100)*a65
a64 += (-459//200)*a65
k1 = integrator.fsalfirst
k2 = f(uprev+dt*a21*k1, p, t+c2*dt)
k3 = f(uprev+dt*(a31*k1+a32*k2), p, t+c3*dt)
k4 = f(uprev+dt*(a41*k1+a42*k2+k3), p, t+c4*dt)
k4 = f(uprev+dt*(a41*k1+a42*k2+2*k3), p, t+c4*dt)
k5 = f(uprev+dt*(a51*k1+a52*k2+a53*k3+a54*k4), p, t+c5*dt)
k6 = f(uprev+dt*(a61*k1+a62*k2+a63*k3+a64*k4+a65*k5), p, t+c6*dt)
u = uprev+dt*(b1*k1+b3*k3+b4*k4+b5*k5+b6*k6)
Expand Down Expand Up @@ -447,7 +450,10 @@ end
w = integrator.alg.w
v = w*dt
## Formula by Z.A. Anastassi, see the Anas5 caches in tableaus/low_order_rk_tableaus.jl for the full citation.
a65 = (-8000//1071)*(-a43*v^5 + 6*tan(v)*v^4 + 24*v^3 - 72*tan(v)*v^2 - 144*v + 144*tan(v))/(v^5*(a43*tan(v)*v + 12 - 10*a43))
a65 = (-8000//1071)*(-a43*(v^5) + 6*tan(v)*(v^4) + 24*(v^3) - 72*tan(v)*(v^2) - 144*v + 144*tan(v))/((v^5)*(a43*tan(v)*v + 12 - 10*a43))
a61 += (-119//200)*a65
a63 += (189//100)*a65
a64 += (-459//200)*a65
@tight_loop_macros for i in uidx
@inbounds tmp[i] = uprev[i]+a21*k1[i]
end
Expand All @@ -457,7 +463,7 @@ end
end
f(k3, tmp, p, t+c3*dt)
@tight_loop_macros for i in uidx
@inbounds tmp[i] = uprev[i]+dt*(a41*k1[i]+a42*k2[i]+k3[i])
@inbounds tmp[i] = uprev[i]+dt*(a41*k1[i]+a42*k2[i]+2*k3[i])
end
f(k4, tmp, p, t+c4*dt)
@tight_loop_macros for i in uidx
Expand Down
16 changes: 8 additions & 8 deletions src/tableaus/low_order_rk_tableaus.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1431,11 +1431,11 @@ Base.@pure function Anas5ConstantCache{T<:CompiledFloats,T2<:CompiledFloats}(::T
a52 = T(2.1375)
a53 = T(-0.2295)
a54 = T(0.401625)
a61 = T(-5.19)
a61 = T(-4.0)
a62 = T(5.0)
a63 = T(1.89)
a64 = T(-2.295)
a65 = T(1)
a63 = T(0.0)
a64 = T(0.0)
a65 = T(0.0)
c2 = T2(0.1)
c3 = T2(0.3333333333333333)
c4 = T2(0.6666666666666667)
Expand All @@ -1460,11 +1460,11 @@ Base.@pure function Anas5ConstantCache{T,T2}(::Type{T},::Type{T2})
a52 = T(171//80)
a53 = T(-459//2000)
a54 = T(3213//8000)
a61 = T(-519//100)
a61 = T(-4//1)
a62 = T(5//1)
a63 = T(189//100)
a64 = T(-459//200)
a65 = T(1//1)
a63 = T(0//1)
a64 = T(0//1)
a65 = T(0//1)
c2 = T2(1//10)
c3 = T2(1//3)
c4 = T2(2//3)
Expand Down
4 changes: 2 additions & 2 deletions test/ode/ode_convergence_tests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -51,8 +51,8 @@ for i = 1:2
@test abs(sim105.𝒪est[:l2]-4) < testTol
sim106 = test_convergence(dts,prob,VCABM5())
@test abs(sim106.𝒪est[:l2]-5) < testTol
sim160 = test_convergence(dts,prob,Anas5(w=5))
@test abs(sim160.𝒪est[:final]-5) < testTol
sim160 = test_convergence(dts,prob,Anas5(w=2))
@test abs(sim160.𝒪est[:l2]-4) < testTol

println("Stiff Solvers")

Expand Down

0 comments on commit 12fcc1b

Please sign in to comment.