Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Fifth order EPIRK methods #417

Merged
merged 10 commits into from
Jul 9, 2018
Merged

Fifth order EPIRK methods #417

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
770d78f
`EPIRK5s3` (out-of-place)
MSeeker1340 Jul 5, 2018
d8012dd
`EPIRK5s3` (in-place)
MSeeker1340 Jul 5, 2018
5669bfe
`EXPRB53s3`
MSeeker1340 Jul 5, 2018
e48544b
Update test script
MSeeker1340 Jul 5, 2018
f4e2884
Merge remote-tracking branch 'origin/master' into expRK
MSeeker1340 Jul 8, 2018
5091d52
`EPIRK5P1` (out-of-place)
MSeeker1340 Jul 8, 2018
dcb53c0
`EPIRK5P1` (in-place)
MSeeker1340 Jul 8, 2018
4263cc7
`EPIRK5P2` (out-of-place)
MSeeker1340 Jul 8, 2018
6c2ed45
`EPIRK5P2` (in-place)
MSeeker1340 Jul 9, 2018
173fefc
Fix `EPIRK5P1` coefficients
MSeeker1340 Jul 9, 2018
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
3 changes: 2 additions & 1 deletion src/OrdinaryDiffEq.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -164,7 +164,8 @@ module OrdinaryDiffEq

export GenericIIF1, GenericIIF2

export LawsonEuler, NorsettEuler, ETD1, ETDRK2, ETDRK3, ETDRK4, HochOst4, Exp4, EPIRK4s3A, EPIRK4s3B, ETD2
export LawsonEuler, NorsettEuler, ETD1, ETDRK2, ETDRK3, ETDRK4, HochOst4, Exp4, EPIRK4s3A, EPIRK4s3B,
EPIRK5s3, EXPRB53s3, EPIRK5P1, EPIRK5P2, ETD2

export SymplecticEuler, VelocityVerlet, VerletLeapfrog, PseudoVerletLeapfrog,
McAte2, Ruth3, McAte3, CandyRoz4, McAte4, McAte42, McAte5,
Expand Down
4 changes: 4 additions & 0 deletions src/alg_utils.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -112,6 +112,10 @@ alg_order(alg::HochOst4) = 4
alg_order(alg::Exp4) = 4
alg_order(alg::EPIRK4s3A) = 4
alg_order(alg::EPIRK4s3B) = 4
alg_order(alg::EPIRK5s3) = 5
alg_order(alg::EPIRK5P1) = 5
alg_order(alg::EPIRK5P2) = 5
alg_order(alg::EXPRB53s3) = 5
alg_order(alg::SplitEuler) = 1
alg_order(alg::ETD2) = 2

Expand Down
2 changes: 1 addition & 1 deletion src/algorithms.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -748,7 +748,7 @@ for Alg in [:LawsonEuler, :NorsettEuler, :ETDRK2, :ETDRK3, :ETDRK4, :HochOst4]
@eval Base.@pure $Alg(;krylov=false, m=30, iop=0) = $Alg(krylov, m, iop)
end
ETD1 = NorsettEuler # alias
for Alg in [:Exp4, :EPIRK4s3A, :EPIRK4s3B]
for Alg in [:Exp4, :EPIRK4s3A, :EPIRK4s3B, :EPIRK5s3, :EXPRB53s3, :EPIRK5P1, :EPIRK5P2]
@eval struct $Alg <: OrdinaryDiffEqExponentialAlgorithm
m::Int
iop::Int
Expand Down
120 changes: 120 additions & 0 deletions src/caches/linear_nonlinear_caches.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -424,6 +424,126 @@ function alg_cache(alg::EPIRK4s3B,u,rate_prototype,uEltypeNoUnits,uBottomEltypeN
EPIRK4s3BCache(u,uprev,tmp,rtmp,rtmp2,K,A,B,KsCache)
end

struct EPIRK5s3ConstantCache <: ExpRKConstantCache end
alg_cache(alg::EPIRK5s3,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,
uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{false}}) = EPIRK5s3ConstantCache()

struct EPIRK5s3Cache{uType,rateType,matType,KsType} <: ExpRKCache
u::uType
uprev::uType
tmp::uType
k::uType
rtmp::rateType
rtmp2::rateType
A::matType
B::matType
KsCache::KsType
end
function alg_cache(alg::EPIRK5s3,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,
tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{true}})
tmp, k = (similar(u) for i = 1:2) # uType caches
rtmp, rtmp2 = (zeros(rate_prototype) for i = 1:2) # rateType caches
YingboMa marked this conversation as resolved.
Show resolved Hide resolved
# Allocate matrices
n = length(u); T = eltype(u)
A = Matrix{T}(n, n) # TODO: sparse Jacobian support
B = zeros(T, n, 5)
YingboMa marked this conversation as resolved.
Show resolved Hide resolved
# Allocate caches for phiv_timestep
maxiter = min(alg.m, n)
KsCache = _phiv_timestep_caches(u, maxiter, 4)
EPIRK5s3Cache(u,uprev,tmp,k,rtmp,rtmp2,A,B,KsCache)
end

struct EXPRB53s3ConstantCache <: ExpRKConstantCache end
alg_cache(alg::EXPRB53s3,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,
uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{false}}) = EXPRB53s3ConstantCache()

struct EXPRB53s3Cache{uType,rateType,matType,KsType} <: ExpRKCache
u::uType
uprev::uType
tmp::uType
rtmp::rateType
rtmp2::rateType
K::matType
A::matType
B::matType
KsCache::KsType
end
function alg_cache(alg::EXPRB53s3,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,
tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{true}})
tmp = similar(u) # uType caches
rtmp, rtmp2 = (zeros(rate_prototype) for i = 1:2) # rateType caches
# Allocate matrices
n = length(u); T = eltype(u)
K = Matrix{T}(n, 2)
A = Matrix{T}(n, n) # TODO: sparse Jacobian support
B = zeros(T, n, 5)
# Allocate caches for phiv_timestep
maxiter = min(alg.m, n)
KsCache = _phiv_timestep_caches(u, maxiter, 4)
EXPRB53s3Cache(u,uprev,tmp,rtmp,rtmp2,K,A,B,KsCache)
end

struct EPIRK5P1ConstantCache <: ExpRKConstantCache end
alg_cache(alg::EPIRK5P1,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,
uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{false}}) = EPIRK5P1ConstantCache()

struct EPIRK5P1Cache{uType,rateType,matType,KsType} <: ExpRKCache
u::uType
uprev::uType
tmp::uType
rtmp::rateType
rtmp2::rateType
K::matType
A::matType
B::matType
KsCache::KsType
end
function alg_cache(alg::EPIRK5P1,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,
tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{true}})
tmp = similar(u) # uType caches
rtmp, rtmp2 = (zeros(rate_prototype) for i = 1:2) # rateType caches
# Allocate matrices
n = length(u); T = eltype(u)
K = Matrix{T}(n, 3)
A = Matrix{T}(n, n) # TODO: sparse Jacobian support
B = zeros(T, n, 4)
# Allocate caches for phiv_timestep
maxiter = min(alg.m, n)
KsCache = _phiv_timestep_caches(u, maxiter, 3)
EPIRK5P1Cache(u,uprev,tmp,rtmp,rtmp2,K,A,B,KsCache)
end

struct EPIRK5P2ConstantCache <: ExpRKConstantCache end
alg_cache(alg::EPIRK5P2,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,
uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{false}}) = EPIRK5P2ConstantCache()

struct EPIRK5P2Cache{uType,rateType,matType,KsType} <: ExpRKCache
u::uType
uprev::uType
tmp::uType
rtmp::rateType
rtmp2::rateType
dR::rateType
K::matType
A::matType
B::matType
KsCache::KsType
end
function alg_cache(alg::EPIRK5P2,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,
tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Type{Val{true}})
tmp = similar(u) # uType caches
rtmp, rtmp2, dR = (zeros(rate_prototype) for i = 1:3) # rateType caches
# Allocate matrices
n = length(u); T = eltype(u)
K = Matrix{T}(n, 3)
A = Matrix{T}(n, n) # TODO: sparse Jacobian support
B = zeros(T, n, 4)
# Allocate caches for phiv_timestep
maxiter = min(alg.m, n)
KsCache = _phiv_timestep_caches(u, maxiter, 3)
EPIRK5P2Cache(u,uprev,tmp,rtmp,rtmp2,dR,K,A,B,KsCache)
end

####################################
# Multistep exponential method caches

Expand Down