Skip to content

Commit

Permalink
Merge pull request #22 from QuEraComputing/john/add-docstrings
Browse files Browse the repository at this point in the history
Better README and added docstrings
  • Loading branch information
johnzl-777 authored Jan 12, 2024
2 parents acb4213 + 0cf1e4d commit 965f5cc
Show file tree
Hide file tree
Showing 3 changed files with 313 additions and 3 deletions.
62 changes: 61 additions & 1 deletion README.md
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

[![codecov](https://codecov.io/github/QuEraComputing/DormandPrince.jl/graph/badge.svg?token=qYZ4V7m0JY)](https://codecov.io/github/QuEraComputing/DormandPrince.jl)

Julia-native implementation of the Dormand-Prince 5th order solver
Julia-native implementation of the Dormand-Prince 5th and 8th order solvers

## Installation

Expand All @@ -22,6 +22,66 @@ DormandPrince is a  
pkg> add DormandPrince
```

## Usage

### Single End-Time

```julia
julia> using DormandPrince

# define your ODE, in this case, dy/dx = 0.85y
julia> function fcn(x, y, f)
f[1] = 0.85y[1]
end

# Create a solver object. We will use the 5th order solver and
# start integrating at t = 0.0 with initial value of 19.0
julia> solver = DP5Solver(fcn, 0.0, [19.0])

# Begin integration up to t = 2.1
julia> integrate!(solver, 2.1)

# get_current_state will return the answer
julia> get_current_state(solver)
```

Both the `DP5Solver` and `DP8Solver`'s are stateful allowing for memory-efficient integration to future end times from the last integrated end point (e.g. if you chose t = 1.0 as your endpoint you can call `integrate!` again with t=2.0 and it will "carry forward" the work starting from t = 1.0 instead of requiring you to set things up all over again).

### Multiple End-Times

```julia
julia> using DormandPrince

# Define ODE
julia> function fcn(x, y, f)
f[1] = 0.85y[1]
end

# Define times of interest to analyze/perform actions on the solution
julia> times = [1.0, 1.1, 1.9, 2.4]

# Create the solver object, with integration starting from t = 0.0 and initial value of 19.0
julia> solver = DP5Solver(fcn, 0.0, [19.0])

# Empty array to store intermediate values
julia> intermediate_values = []

# Use callback to store intermediate values. The solver object will also be mutated to store the solution
# found at the last time point.
julia> integrate!(solver, times) do time, val
push!(intermediate_values, val[])
end

```

You may also create a `SolverIterator` that can then be iterated over. Note that this will require a fresh solver

```julia
for (time, value) in SolverIterator(solver, times)
println("Time: ", time, " ", "Value: ", value[])
end
```

## License

Apache License 2.0
81 changes: 81 additions & 0 deletions src/interface.jl
Original file line number Diff line number Diff line change
@@ -1,4 +1,28 @@
"""
struct SolverIterator{T <: Real}
SolverIterator(solver, times)
Given a solver and a vector of times,
this iterator will return the state of the solver at each time.
The solver will be mutated to hold the solution from the last time given.
# Examples
```julia
julia> solver = DP5Solver(fcn, 0.0, [0.0])
julia> times = [1.0, 2.0, 3.0]
julia> intermediate_vals = []
julia> for (time, value) in SolverIterator(solver, times)
push!(intermediate_values, value[])
end
```
"""
struct SolverIterator{T <: Real}
solver::AbstractDPSolver{T}
times::AbstractVector{T}
Expand Down Expand Up @@ -32,12 +56,69 @@ end
get_current_state(::AbstractDPSolver) = error("not implemented")
integrate_core!(::AbstractDPSolver{T}, ::T) where T = error("not implemented")

"""
integrate!(solver, time)
Integrate the ODE problem that is part of solver to the end time `time`.
`solver` can be a `DP5Solver` or a `DP8Solver` type.
The `solver` is mutated to hold the solution and `solver` state at the time `time`, allowing for
efficient integration to future end times of interest through subsequent calls to `integrate!` with
later times.
# Examples
```julia
julia> function fcn(x, y, f)
f[1] = 0.85y[1]
end
julia> solver = DP5Solver(fcn, 0.0, [19.0])
# Integrate from initial time of 0.0 to 1.0
julia> integrate!(solver, 1.0)
# integrate from the last time of 1.0 to 2.0
julia> integrate!(solver, 2.0)
# Solver object now holds solution and state at 2.0
julia> get_current_state(solver)
```
"""
function integrate!(solver::AbstractDPSolver{T}, time::T) where T <: Real
report = integrate_core!(solver, time)
if report.idid != COMPUTATION_SUCCESSFUL
throw(DPException(report))
end
end

"""
integrate!(callback, solver, times)
Integrate the ODE problem that is part of `solver` to the end times `times` and apply the callback on
the time and solution at each time.
`solver` can be a `DP5Solver` or a `DP8Solver` type.
At the end of all the times the solver holds the solution and solver state at the last time in `times`.
# Examples
```julia
julia> function fcn(x, y, f)
f[1] = 0.85y[1]
end
julia> times = [1.0, 1.1, 1.9, 2.4]
julia> solver = DP5Solver(fcn, 0.0, [19.0])
julia> intermediate_values = []
julia> integrate!(solver, times) do time, val
push!(intermediate_values, val[])
end
```
"""
function integrate!(callback, solver::AbstractDPSolver{T}, times::AbstractVector{T}; sort_times::Bool = true) where {T <: Real}
times = sort_times ? sort(collect(times)) : times

Expand Down
Loading

0 comments on commit 965f5cc

Please sign in to comment.