Add defect control adaptivity

.gitignore

@@ -2,3 +2,5 @@
 *.jl.*.cov
 *.jl.mem
 Manifest.toml
+.vscode
+.vscode/*

src/BoundaryValueDiffEq.jl

@@ -8,22 +8,43 @@ import ForwardDiff, BandedMatrices, FiniteDiff
 
 import TruncatedStacktraces: @truncate_stacktrace
 
-struct MIRKTableau{T, cType, vType, bType, xType}
+struct MIRKTableau{T, cType, vType, bType, xType, sType, starType, tauType}
     c::cType
     v::vType
     b::bType
     x::xType
+    """Discrete stages of MIRK formula"""
+    s::sType
+    """Number of stages to form the interpolant"""
+    s_star::starType
+    """Defect sampling point"""
+    tau::tauType
 end
 
-function MIRKTableau(c, v, b, x)
-    @assert eltype(c) == eltype(v) == eltype(b) == eltype(x)
-    return MIRKTableau{eltype(c), typeof(c), typeof(v), typeof(b), typeof(x)}(c, v, b, x)
+function MIRKTableau(c, v, b, x, s, s_star, tau)
+    @assert eltype(c) == eltype(v) == eltype(b) == eltype(x) == eltype(tau)
+    return MIRKTableau{
+        eltype(c),
+        typeof(c),
+        typeof(v),
+        typeof(b),
+        typeof(x),
+        typeof(s),
+        typeof(s_star),
+        typeof(tau),
+    }(c,
+        v,
+        b,
+        x,
+        s,
+        s_star,
+        tau)
 end
 
 @truncate_stacktrace MIRKTableau 1
 
 # ODE BVP problem system
-struct BVPSystem{T, U <: AbstractArray, P, F, B, S}
+mutable struct BVPSystem{T, U <: AbstractArray, P, F, B, S}
     order::Int                  # The order of MIRK method
     M::Int                      # Number of equations in the ODE system
     N::Int                      # Number of nodes in the mesh
@@ -46,6 +67,7 @@ include("cache.jl")
 include("collocation.jl")
 include("jacobian.jl")
 include("solve.jl")
+include("adaptivity.jl")
 
 export Shooting
 export GeneralMIRK4, GeneralMIRK5, GeneralMIRK6