Merge pull request #11 from a-ev/docstrings

Working on docstrings

.github/workflows/documentation.yml

@@ -18,8 +18,6 @@ jobs:
       - uses: julia-actions/setup-julia@v2
         with:
           version: '1.10'
-      - run: julia --project -e 'using Pkg; Pkg.add(url="https://github.com/chrisdjscott/DACE_jll.jl.git")'
-        shell: bash
       - name: Install dependencies
         run: julia --project=docs/ -e 'using Pkg; Pkg.instantiate()'
       - name: Build and deploy

.github/workflows/test.yml

@@ -1,6 +1,6 @@
 name: Run tests
 
-on: 
+on:
   push:
   pull_request:
 
@@ -22,8 +22,6 @@ jobs:
         with:
           version: ${{ matrix.julia-version }}
           arch: ${{ matrix.julia-arch }}
-      - run: julia --project -e 'using Pkg; Pkg.add(url="https://github.com/chrisdjscott/DACE_jll.jl.git")'
-        shell: bash
       - uses: julia-actions/julia-buildpkg@v1
       - uses: julia-actions/julia-runtest@v1
         # with:

Project.toml

@@ -7,10 +7,15 @@ version = "0.1.0"
 CxxWrap = "1f15a43c-97ca-5a2a-ae31-89f07a497df4"
 DACE_jll = "40de70a5-cf80-53d9-bda4-3aa67fea2f4f"
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
+SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
+DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 
 [compat]
-CxxWrap = "0.15"
-DACE_jll = ">= 0.0.11"
+CxxWrap = "~0.16"
+DACE_jll = "~0.4"
+DiffEqBase = "^6"
+SpecialFunctions = "^2"
+DocStringExtensions = "^0.9"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

docs/Project.toml

@@ -1,3 +1,16 @@
 [deps]
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
+OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+CxxWrap = "1f15a43c-97ca-5a2a-ae31-89f07a497df4"
+DACE_jll = "40de70a5-cf80-53d9-bda4-3aa67fea2f4f"
+DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
+SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
+DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
+
+[compat]
+CxxWrap = "~0.16"
+DACE_jll = "~0.4"
+DiffEqBase = "^6"
+SpecialFunctions = "^2"
+DocStringExtensions = "^0.9"

docs/make.jl

@@ -12,6 +12,7 @@ const OUTPUT_DIR = joinpath(@__DIR__, "src", "generated")
 example_scripts = [
     "sine.jl",
     "polynomial_inversion.jl",
+    "ode_integration.jl"
 ]
 
 println("Building examples...")
@@ -29,6 +30,7 @@ println("Finished building examples")
 example_pages = [
     "Sine function" => "generated/sine.md",
     "Polynomial inversion" => "generated/polynomial_inversion.md",
+    "ODE integration" => "generated/ode_integration.md"
 ]
 
 tutorial_pages = [
@@ -42,6 +44,7 @@ pages = [
 #    "Developer Guide" => "developing.md",
     "Tutorials" => tutorial_pages,
     "Examples" => example_pages,
+    "API" => "api.md"
 ]
 
 # build and deploy docs

docs/src/api.md

@@ -0,0 +1,12 @@
+# API
+
+```@meta
+CurrentModule = DACE
+```
+
+```@index
+```
+
+```@autodocs
+Modules = [DACE]
+```

examples/ode_integration.jl

@@ -1,63 +1,72 @@
-"""
-This tutorials shows how to perform the numerical 
-integration of an ordinary differential equation (ODE) 
-using DA objects as the state components' type.
+# # ODE integration example
+#
+# This tutorial shows how to perform the numerical integration of an ordinary differential 
+# equation (ODE) using DA objects as the state components' type.
+# It also demonstrates how to extract the state transition matrix (STM) from the polynomial 
+# expansion of the dynamics flow.
+#
+# In this case, the state vector represents an object in orbit around a central body and 
+# subject only to the gravitational pull of the latter.
+# Its motion is described by the differential equations for the Kepler problem expressed in 
+# Cartesian coordinates.
 
-It propagates the differential equations for 
-the Keplerian motion in Cartesian parameters.
-"""
+# ## Install dependencies
+#
+# Make sure the required packages are installed
+
+# ```julia
+# using Pkg
+# Pkg.add("https://github.com/a-ev/DACE.jl.git")
+# ```
+
+# ## Using DACE
+#
+# Load the required modules
 
 using OrdinaryDiffEq
 using DACE
 
-# initial parameters defining a circular orbit propagated for 1 revolution
+# Define the parameters and intial conditions for a circular orbit with a normalized radius 
+# equal to one
 μ = 1.0
-ts = 0.0
+x0 = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0]
+
+# Set the integration time span equal to one revolution
+t0 = 0.0
 tf = 2π
-xs = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0]
 
-# equations of motion
-function kepler_ode!(du, u, p, t)
+# Define the equations of motion for the resticted two-body problem
+function kepler_ode!(du, u, μ, _)
     du[1:3] .= u[4:6]
     du[4:6] .= -u[1:3] * μ ./ sum(u[1:3].^2)^(3/2)
 end
 
-# nominal solution
-prob = ODEProblem(kepler_ode!, xs, [ts, tf])
-sol = solve(prob, DP8(), abstol=1e-12, reltol=1e-12)
+# Compute the nominal solution
+prob = ODEProblem(kepler_ode!, x0, [t0, tf], μ)
+sol = solve(prob, Vern9(), abstol=1e-12, reltol=1e-12)
+
+xf = sol.u[end]
 
-# DACE initialization for computations at second order
+# Initialize DACE to compute second-order expansions
 DACE.init(2,6)
 
-# initial conditions as DA vector
+# Define the initial conditions as a DA vector
 dx_dace = [DA(i,1.0) for i in 1:6]
-xs_dace = xs .+ dx_dace
+x0_dace = x0 .+ dx_dace
 
-# polynomial expansion of the dynamics flow
-sol_dace = solve(remake(prob, u0=xs_dace), DP8(), abstol=1e-12, reltol=1e-12)
+println(x0_dace)
 
-# verify that the time steps chosen by the integrator coincide
-t_dace = sol_dace.t
+# Compute the polynomial expansion of the dynamics flow
+sol_dace = solve(remake(prob, u0=x0_dace), Vern9(), abstol=1e-12, reltol=1e-12)
 
-println("Time steps coincide: " * string(all(t_dace .≈ sol.t)))
-
-# verify the (nominal) final state
 xf_dace = sol_dace.u[end]
 xf_cons = DACE.cons.(xf_dace)
 
-println("Max abs error on final state: " * string(maximum(abs.(xf_cons - sol.u[end]))))
-
-# extract and print the state transition matrix (STM)
-function dace_jacobian(v::Vector{<:DA})
-    jac = Matrix{Float64}(undef, length(v), DACE.getMaxVariables())
-    for i in 1:length(v)
-        for j in 1:DACE.getMaxVariables()
-            jac[i,j] = DACE.cons(DACE.deriv(v[i],j))
-        end
-    end
-    return jac
-end;
-
-stm_dace = dace_jacobian(xf_dace)
-println("State transition matrix:")
-display(stm_dace)
+# Verify the (nominal) final state
+
+println("Max abs error on final state: " * string(maximum(abs.(xf_cons - xf))))
+
+# Extract the state transition matrix (STM)
+
+stm_dace = DACE.jacobian(xf_dace)
+stm_cons = DACE.cons.(stm_dace)