From 5a21c14eaf39b809024dac2f5a23d5d19cb3e8f2 Mon Sep 17 00:00:00 2001 From: Olivier Cots Date: Thu, 11 Jun 2026 22:48:42 +0200 Subject: [PATCH 1/5] [WIP] Placeholder for issue #445 - integrate utils from CTModels From e801b6b7604861d1b6349361d1f4e2192162f7f2 Mon Sep 17 00:00:00 2001 From: Olivier Cots Date: Thu, 11 Jun 2026 23:03:40 +0200 Subject: [PATCH 2/5] Add Interpolation module and extend Core with utilities from CTModels - New Interpolation module with ctinterpolate functions - Extended Core with matrix_utils, function_utils, macros - Updated docs/api_reference.jl for new modules - Moved DeepWiki tip in index.md --- docs/api_reference.jl | 18 ++++++++++++++++++ docs/src/index.md | 8 ++++---- src/CTBase.jl | 4 ++++ src/Core/Core.jl | 9 ++++++++- 4 files changed, 34 insertions(+), 5 deletions(-) diff --git a/docs/api_reference.jl b/docs/api_reference.jl index 0904cbc7..d8a8d4a6 100644 --- a/docs/api_reference.jl +++ b/docs/api_reference.jl @@ -22,6 +22,9 @@ function generate_api_reference(src_dir::String) joinpath("Core", "Core.jl"), joinpath("Core", "default.jl"), joinpath("Core", "types.jl"), + joinpath("Core", "matrix_utils.jl"), + joinpath("Core", "function_utils.jl"), + joinpath("Core", "macros.jl"), ), ], exclude=EXCLUDE_SYMBOLS, @@ -31,6 +34,21 @@ function generate_api_reference(src_dir::String) title_in_menu="Core", filename="core", ), + CTBase.automatic_reference_documentation(; + subdirectory="api", + primary_modules=[ + CTBase.Interpolation => src( + joinpath("Interpolation", "Interpolation.jl"), + joinpath("Interpolation", "ctinterpolate.jl"), + ), + ], + exclude=EXCLUDE_SYMBOLS, + public=true, + private=false, + title="Interpolation", + title_in_menu="Interpolation", + filename="interpolation", + ), CTBase.automatic_reference_documentation(; subdirectory="api", primary_modules=[ diff --git a/docs/src/index.md b/docs/src/index.md index a92bc260..d9e5b1c2 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -1,9 +1,5 @@ # CTBase.jl — Ecosystem Foundation -!!! tip "Ask DeepWiki" - [![Ask DeepWiki](https://deepwiki.com/badge.svg)](https://deepwiki.com/control-toolbox/CTBase.jl) offers an interactive, AI-generated overview of this codebase. Answers may be inaccurate — use this - reference documentation as the source of truth. - ```@meta CurrentModule = CTBase ``` @@ -20,6 +16,10 @@ It provides the **base layer** shared by all packages: common types, structured Downstream packages (e.g. [OptimalControl.jl](https://github.com/control-toolbox/OptimalControl.jl)) may re-export selected symbols for convenience. +!!! tip "Ask DeepWiki" + [![Ask DeepWiki](https://deepwiki.com/badge.svg)](https://deepwiki.com/control-toolbox/CTBase.jl) offers an interactive, AI-generated overview of this codebase. Answers may be inaccurate — use this + reference documentation as the source of truth. + ## Submodule overview | Submodule | Role | diff --git a/src/CTBase.jl b/src/CTBase.jl index e1ddd1b7..46466231 100644 --- a/src/CTBase.jl +++ b/src/CTBase.jl @@ -32,4 +32,8 @@ using .Descriptions include(joinpath(@__DIR__, "Extensions", "Extensions.jl")) using .Extensions +# Interpolation module - interpolation utilities +include(joinpath(@__DIR__, "Interpolation", "Interpolation.jl")) +using .Interpolation + end diff --git a/src/Core/Core.jl b/src/Core/Core.jl index 2a6f0454..75e37d97 100644 --- a/src/Core/Core.jl +++ b/src/Core/Core.jl @@ -13,7 +13,14 @@ import DocStringExtensions: TYPEDEF, TYPEDSIGNATURES include("types.jl") include("default.jl") +# Private utilities +include("function_utils.jl") +include("macros.jl") + +# Public utilities +include("matrix_utils.jl") + # Export public API -export ctNumber +export ctNumber, matrix2vec, to_out_of_place, @ensure end # module From cebb4e394525cc0ddf042f3b3b842441f26c7ef5 Mon Sep 17 00:00:00 2001 From: Olivier Cots Date: Thu, 11 Jun 2026 23:11:06 +0200 Subject: [PATCH 3/5] Add Interpolation module and Core utilities from CTModels - Interpolation module: ctinterpolate, ctinterpolate_constant - Core utilities: matrix_utils, function_utils, macros - Tests for all new modules --- src/Core/function_utils.jl | 32 +++++ src/Core/macros.jl | 26 ++++ src/Core/matrix_utils.jl | 39 ++++++ src/Interpolation/Interpolation.jl | 18 +++ src/Interpolation/ctinterpolate.jl | 87 +++++++++++++ test/suite/core/test_function_utils.jl | 74 +++++++++++ test/suite/core/test_macros.jl | 64 +++++++++ test/suite/core/test_matrix_utils.jl | 72 +++++++++++ .../suite/interpolation/test_interpolation.jl | 122 ++++++++++++++++++ 9 files changed, 534 insertions(+) create mode 100644 src/Core/function_utils.jl create mode 100644 src/Core/macros.jl create mode 100644 src/Core/matrix_utils.jl create mode 100644 src/Interpolation/Interpolation.jl create mode 100644 src/Interpolation/ctinterpolate.jl create mode 100644 test/suite/core/test_function_utils.jl create mode 100644 test/suite/core/test_macros.jl create mode 100644 test/suite/core/test_matrix_utils.jl create mode 100644 test/suite/interpolation/test_interpolation.jl diff --git a/src/Core/function_utils.jl b/src/Core/function_utils.jl new file mode 100644 index 00000000..effdfe26 --- /dev/null +++ b/src/Core/function_utils.jl @@ -0,0 +1,32 @@ +""" +$(TYPEDSIGNATURES) + +Convert an in-place function `f!` to an out-of-place function `f`. + +The resulting function `f` returns a vector of type `T` and length `n` by first allocating +memory and then calling `f!` to fill it. + +# Arguments +- `f!`: An in-place function of the form `f!(result, args...)`. +- `n`: The length of the output vector. +- `T`: The element type of the output vector (default is `Float64`). + +# Returns +An out-of-place function `f(args...; kwargs...)` that returns the result as a vector or +scalar, depending on `n`. + +# Example +```julia-repl +julia> f!(r, x) = (r[1] = sin(x); r[2] = cos(x)) +julia> f = to_out_of_place(f!, 2) +julia> f(π/4) # returns approximately [0.707, 0.707] +``` +""" +function to_out_of_place(f!, n; T=Float64) + function f(args...; kwargs...) + r = zeros(T, n) + f!(r, args...; kwargs...) + return n == 1 ? r[1] : r + end + return isnothing(f!) ? nothing : f +end diff --git a/src/Core/macros.jl b/src/Core/macros.jl new file mode 100644 index 00000000..7dc6fc72 --- /dev/null +++ b/src/Core/macros.jl @@ -0,0 +1,26 @@ +""" + @ensure condition exception + +Throws the provided `exception` if `condition` is false. + +# Usage +```julia-repl +julia> @ensure true Exceptions.IncorrectArgument("This won't throw") +julia> @ensure false Exceptions.IncorrectArgument("This will throw") +ERROR: IncorrectArgument("This will throw") +``` + +# Arguments +- `condition`: A Boolean expression to test. +- `exception`: An instance of an exception to throw if `condition` is false. + +# Throws +- The provided `exception` if the condition is not satisfied. +""" +macro ensure(cond, exc) + return esc(:( + if !($cond) + throw($exc) + end + )) +end diff --git a/src/Core/matrix_utils.jl b/src/Core/matrix_utils.jl new file mode 100644 index 00000000..553ad373 --- /dev/null +++ b/src/Core/matrix_utils.jl @@ -0,0 +1,39 @@ +""" +$(TYPEDSIGNATURES) + +Return the default dimension for matrix storage. +""" +__matrix_dimension_storage() = 1 + +""" +$(TYPEDSIGNATURES) + +Transform a matrix into a vector of vectors along the specified dimension. + +Each row or column of the matrix `A` is extracted and stored as an individual vector, +depending on `dim`. + +# Arguments +- `A`: A matrix of elements of type `<:ctNumber`. +- `dim`: The dimension along which to split the matrix (`1` for rows, `2` for columns). + Defaults to `1`. + +# Returns +A `Vector` of `Vector`s extracted from the rows or columns of `A`. + +# Note +This is useful when data needs to be represented as a sequence of state or control vectors +in optimal control problems. + +# Example +```julia-repl +julia> A = [1 2 3; 4 5 6] +julia> matrix2vec(A, 1) # splits into rows: [[1, 2, 3], [4, 5, 6]] +julia> matrix2vec(A, 2) # splits into columns: [[1, 4], [2, 5], [3, 6]] +``` +""" +function matrix2vec( + A::Matrix{<:ctNumber}, dim::Int=__matrix_dimension_storage() +)::Vector{<:Vector{<:ctNumber}} + return dim == 1 ? [A[i, :] for i in 1:size(A, 1)] : [A[:, i] for i in 1:size(A, 2)] +end diff --git a/src/Interpolation/Interpolation.jl b/src/Interpolation/Interpolation.jl new file mode 100644 index 00000000..03d0fc79 --- /dev/null +++ b/src/Interpolation/Interpolation.jl @@ -0,0 +1,18 @@ +""" + Interpolation + +Interpolation utilities for the Control Toolbox (CT) ecosystem. + +# Public API +- `ctinterpolate`: linear interpolation with flat extrapolation +- `ctinterpolate_constant`: piecewise-constant (steppost) interpolation +""" +module Interpolation + +import DocStringExtensions: TYPEDSIGNATURES + +include(joinpath(@__DIR__, "ctinterpolate.jl")) + +export ctinterpolate, ctinterpolate_constant + +end # module Interpolation diff --git a/src/Interpolation/ctinterpolate.jl b/src/Interpolation/ctinterpolate.jl new file mode 100644 index 00000000..35c4a749 --- /dev/null +++ b/src/Interpolation/ctinterpolate.jl @@ -0,0 +1,87 @@ +""" +$(TYPEDSIGNATURES) + +Return a linear interpolation function for the data `f` defined at points `x`. + +This function creates a one-dimensional linear interpolant with flat extrapolation +beyond the bounds of `x` (returns `f[1]` for `t < x[1]` and `f[end]` for `t >= x[end]`). + +# Arguments +- `x`: A vector of points at which the values `f` are defined. +- `f`: A vector of values to interpolate. + +# Returns +A callable interpolation object that can be evaluated at new points. + +# Example +```julia-repl +julia> x = [0.0, 1.0, 2.0] +julia> f = [1.0, 2.0, 3.0] +julia> interp = ctinterpolate(x, f) +julia> interp(0.5) # Returns 1.5 (linear interpolation) +julia> interp(-1.0) # Returns 1.0 (flat extrapolation) +julia> interp(3.0) # Returns 3.0 (flat extrapolation) +``` +""" +function ctinterpolate(x, f) + function linear_interp(t) + if t < x[1] + return f[1] + elseif t >= x[end] + return f[end] + else + i = searchsortedlast(x, t) + if i == length(x) + return f[end] + end + α = (t - x[i]) / (x[i + 1] - x[i]) + return f[i] + α * (f[i + 1] - f[i]) + end + end + return linear_interp +end + +""" +$(TYPEDSIGNATURES) + +Return a piecewise-constant interpolation function for the data `f` defined at points `x`. + +This function creates a right-continuous piecewise constant interpolant: +the value at knot `x[i]` is held constant on the interval `[x[i], x[i+1})`. + +This implements the standard steppost behavior for optimal control: +- `u(t_i) = u_i` (value at the knot) +- `u(t) = u_i` for all `t ∈ [t_i, t_{i+1})` +- Right-continuous: `lim_{t→t_i^+} u(t) = u(t_i)` + +# Arguments +- `x`: A vector of points at which the values `f` are defined. +- `f`: A vector of values to interpolate. + +# Returns +A callable interpolation object that can be evaluated at new points. + +# Example +```julia-repl +julia> x = [0.0, 1.0, 2.0] +julia> f = [1.0, 2.0, 3.0] +julia> interp = ctinterpolate_constant(x, f) +julia> interp(0.0) # Returns 1.0 (value at x[1]) +julia> interp(0.5) # Returns 1.0 (held from x[1] on [0.0, 1.0)) +julia> interp(1.0) # Returns 2.0 (value at x[2], right-continuous) +julia> interp(1.5) # Returns 2.0 (held from x[2] on [1.0, 2.0)) +``` +""" +function ctinterpolate_constant(x, f) + function steppost_interp(t) + if t < x[1] + return f[1] + elseif t >= x[end] + return f[end] + else + i = searchsortedlast(x, t) + return f[i] + end + end + return steppost_interp +end diff --git a/test/suite/core/test_function_utils.jl b/test/suite/core/test_function_utils.jl new file mode 100644 index 00000000..3f82837d --- /dev/null +++ b/test/suite/core/test_function_utils.jl @@ -0,0 +1,74 @@ +module TestCoreFunctionUtils + +import Test +import CTBase.Core + +const VERBOSE = isdefined(Main, :TestOptions) ? Main.TestOptions.VERBOSE : true +const SHOWTIMING = isdefined(Main, :TestOptions) ? Main.TestOptions.SHOWTIMING : true + +function test_function_utils() + Test.@testset verbose = VERBOSE showtiming = SHOWTIMING "to_out_of_place" begin + + Test.@testset "to_out_of_place - basic conversion" begin + f!(r, x) = (r[1] = sin(x); r[2] = cos(x)) + f = Core.to_out_of_place(f!, 2) + result = f(π / 4) + Test.@test result isa Vector + Test.@test length(result) == 2 + Test.@test result[1] ≈ sin(π / 4) + Test.@test result[2] ≈ cos(π / 4) + end + + Test.@testset "to_out_of_place - scalar output (n=1)" begin + g!(r, x) = (r[1] = x^2) + g = Core.to_out_of_place(g!, 1) + result = g(3.0) + Test.@test result isa Float64 + Test.@test result ≈ 9.0 + end + + Test.@testset "to_out_of_place - with kwargs" begin + h!(r, x; scale=1.0) = (r[1] = x * scale; r[2] = x^2 * scale) + h = Core.to_out_of_place(h!, 2) + Test.@test h(2.0)[1] ≈ 2.0 + Test.@test h(2.0)[2] ≈ 4.0 + Test.@test h(2.0; scale=3.0)[1] ≈ 6.0 + Test.@test h(2.0; scale=3.0)[2] ≈ 12.0 + end + + Test.@testset "to_out_of_place - multiple arguments" begin + k!(r, x, y) = (r[1] = x + y; r[2] = x * y) + k = Core.to_out_of_place(k!, 2) + result = k(3.0, 4.0) + Test.@test result[1] ≈ 7.0 + Test.@test result[2] ≈ 12.0 + end + + Test.@testset "to_out_of_place - custom type" begin + m!(r, x) = (r[1] = x + 1; r[2] = x + 2) + m = Core.to_out_of_place(m!, 2; T=Int) + result = m(5) + Test.@test result isa Vector{Int} + Test.@test result[1] == 6 + Test.@test result[2] == 7 + end + + Test.@testset "to_out_of_place - nothing input" begin + Test.@test Core.to_out_of_place(nothing, 2) === nothing + end + + Test.@testset "to_out_of_place - larger output" begin + big!(r, x) = (for i in 1:5; r[i] = x * i; end) + big = Core.to_out_of_place(big!, 5) + result = big(2.0) + Test.@test length(result) == 5 + Test.@test result == [2.0, 4.0, 6.0, 8.0, 10.0] + end + + end + return nothing +end + +end # module TestCoreFunctionUtils + +test_function_utils() = TestCoreFunctionUtils.test_function_utils() diff --git a/test/suite/core/test_macros.jl b/test/suite/core/test_macros.jl new file mode 100644 index 00000000..839890ad --- /dev/null +++ b/test/suite/core/test_macros.jl @@ -0,0 +1,64 @@ +module TestCoreMacros + +import Test +import CTBase.Core +import CTBase.Exceptions + +const VERBOSE = isdefined(Main, :TestOptions) ? Main.TestOptions.VERBOSE : true +const SHOWTIMING = isdefined(Main, :TestOptions) ? Main.TestOptions.SHOWTIMING : true + +function test_macros() + Test.@testset verbose = VERBOSE showtiming = SHOWTIMING "@ensure" begin + + Test.@testset "@ensure - condition true" begin + x = 5 + Test.@test_nowarn Core.@ensure x > 0 Exceptions.IncorrectArgument("x must be positive") + Test.@test_nowarn Core.@ensure x == 5 Exceptions.IncorrectArgument("x must be 5") + end + + Test.@testset "@ensure - condition false" begin + x = -5 + Test.@test_throws Exceptions.IncorrectArgument Core.@ensure x > 0 Exceptions.IncorrectArgument("x must be positive") + y = 10 + Test.@test_throws Exceptions.IncorrectArgument Core.@ensure y < 0 Exceptions.IncorrectArgument("y must be negative") + end + + Test.@testset "@ensure - with different exception types" begin + x = 0 + Test.@test_throws ArgumentError Core.@ensure x != 0 ArgumentError("x cannot be zero") + Test.@test_throws DomainError Core.@ensure x > 0 DomainError(x, "x must be positive") + end + + Test.@testset "@ensure - complex conditions" begin + x = 5 + y = 10 + Test.@test_nowarn Core.@ensure x < y Exceptions.IncorrectArgument("x must be less than y") + Test.@test_throws Exceptions.IncorrectArgument Core.@ensure x > y Exceptions.IncorrectArgument("x must be greater than y") + Test.@test_nowarn Core.@ensure (x > 0 && y > 0) Exceptions.IncorrectArgument("both must be positive") + Test.@test_throws Exceptions.IncorrectArgument Core.@ensure (x < 0 || y < 0) Exceptions.IncorrectArgument("at least one must be negative") + end + + Test.@testset "@ensure - exception message verification" begin + x = -5 + try + Core.@ensure x > 0 Exceptions.IncorrectArgument("x must be positive") + Test.@test false + catch e + Test.@test e isa Exceptions.IncorrectArgument + Test.@test e.msg == "x must be positive" + end + end + + Test.@testset "@ensure - with type checks" begin + x = 5 + Test.@test_nowarn Core.@ensure x isa Int Exceptions.IncorrectArgument("x must be an Int") + Test.@test_throws Exceptions.IncorrectArgument Core.@ensure x isa String Exceptions.IncorrectArgument("x must be a String") + end + + end + return nothing +end + +end # module TestCoreMacros + +test_macros() = TestCoreMacros.test_macros() diff --git a/test/suite/core/test_matrix_utils.jl b/test/suite/core/test_matrix_utils.jl new file mode 100644 index 00000000..a517b03a --- /dev/null +++ b/test/suite/core/test_matrix_utils.jl @@ -0,0 +1,72 @@ +module TestCoreMatrixUtils + +import Test +import CTBase.Core + +const VERBOSE = isdefined(Main, :TestOptions) ? Main.TestOptions.VERBOSE : true +const SHOWTIMING = isdefined(Main, :TestOptions) ? Main.TestOptions.SHOWTIMING : true + +function test_matrix_utils() + Test.@testset verbose = VERBOSE showtiming = SHOWTIMING "matrix2vec" begin + + Test.@testset "matrix2vec - dimension 1 (rows)" begin + A = [0 1; 2 3] + V = Core.matrix2vec(A) + Test.@test V isa Vector{<:Vector} + Test.@test length(V) == 2 + Test.@test V[1] == [0, 1] + Test.@test V[2] == [2, 3] + V1 = Core.matrix2vec(A, 1) + Test.@test V1 == V + end + + Test.@testset "matrix2vec - dimension 2 (columns)" begin + A = [0 1; 2 3] + W = Core.matrix2vec(A, 2) + Test.@test W isa Vector{<:Vector} + Test.@test length(W) == 2 + Test.@test W[1] == [0, 2] + Test.@test W[2] == [1, 3] + end + + Test.@testset "matrix2vec - larger matrix" begin + B = [1 2 3; 4 5 6] + rows = Core.matrix2vec(B, 1) + Test.@test length(rows) == 2 + Test.@test rows[1] == [1, 2, 3] + Test.@test rows[2] == [4, 5, 6] + cols = Core.matrix2vec(B, 2) + Test.@test length(cols) == 3 + Test.@test cols[1] == [1, 4] + Test.@test cols[2] == [2, 5] + Test.@test cols[3] == [3, 6] + end + + Test.@testset "matrix2vec - single row/column" begin + R = [1 2 3] + Test.@test length(Core.matrix2vec(R, 1)) == 1 + Test.@test Core.matrix2vec(R, 1)[1] == [1, 2, 3] + Test.@test length(Core.matrix2vec(R, 2)) == 3 + C = reshape([1, 2, 3], 3, 1) + Test.@test length(Core.matrix2vec(C, 1)) == 3 + Test.@test length(Core.matrix2vec(C, 2)) == 1 + Test.@test Core.matrix2vec(C, 2)[1] == [1, 2, 3] + end + + Test.@testset "matrix2vec - Float64 matrix" begin + F = [1.5 2.5; 3.5 4.5] + V = Core.matrix2vec(F, 1) + Test.@test V[1] ≈ [1.5, 2.5] + Test.@test V[2] ≈ [3.5, 4.5] + W = Core.matrix2vec(F, 2) + Test.@test W[1] ≈ [1.5, 3.5] + Test.@test W[2] ≈ [2.5, 4.5] + end + + end + return nothing +end + +end # module TestCoreMatrixUtils + +test_matrix_utils() = TestCoreMatrixUtils.test_matrix_utils() diff --git a/test/suite/interpolation/test_interpolation.jl b/test/suite/interpolation/test_interpolation.jl new file mode 100644 index 00000000..96544431 --- /dev/null +++ b/test/suite/interpolation/test_interpolation.jl @@ -0,0 +1,122 @@ +module TestInterpolation + +import Test +import CTBase.Interpolation + +const VERBOSE = isdefined(Main, :TestOptions) ? Main.TestOptions.VERBOSE : true +const SHOWTIMING = isdefined(Main, :TestOptions) ? Main.TestOptions.SHOWTIMING : true + +function test_interpolation() + Test.@testset verbose = VERBOSE showtiming = SHOWTIMING "Interpolation" begin + + Test.@testset "ctinterpolate - basic linear interpolation" begin + x = [0.0, 1.0, 2.0] + f = [0.0, 1.0, 0.0] + interp = Interpolation.ctinterpolate(x, f) + Test.@test interp(0.0) ≈ 0.0 + Test.@test interp(1.0) ≈ 1.0 + Test.@test interp(2.0) ≈ 0.0 + Test.@test interp(0.5) ≈ 0.5 + Test.@test interp(1.5) ≈ 0.5 + end + + Test.@testset "ctinterpolate - extrapolation" begin + x = [0.0, 1.0, 2.0] + f = [1.0, 2.0, 3.0] + interp = Interpolation.ctinterpolate(x, f) + Test.@test interp(-1.0) ≈ 1.0 + Test.@test interp(3.0) ≈ 3.0 + end + + Test.@testset "ctinterpolate - constant function" begin + x = [0.0, 1.0, 2.0, 3.0] + f = [5.0, 5.0, 5.0, 5.0] + interp = Interpolation.ctinterpolate(x, f) + Test.@test interp(0.5) ≈ 5.0 + Test.@test interp(1.5) ≈ 5.0 + Test.@test interp(2.5) ≈ 5.0 + end + + Test.@testset "ctinterpolate - non-uniform grid" begin + x = [0.0, 0.1, 0.5, 1.0, 2.0] + f = [0.0, 1.0, 2.0, 3.0, 4.0] + interp = Interpolation.ctinterpolate(x, f) + Test.@test interp(0.05) ≈ 0.5 + Test.@test interp(0.3) ≈ 1.5 + Test.@test interp(1.5) ≈ 3.5 + end + + Test.@testset "ctinterpolate - vector values" begin + x = [0.0, 1.0, 2.0] + f = [[0.0, 0.0], [1.0, 2.0], [2.0, 4.0]] + interp = Interpolation.ctinterpolate(x, f) + Test.@test interp(0.0) ≈ [0.0, 0.0] + Test.@test interp(1.0) ≈ [1.0, 2.0] + Test.@test interp(2.0) ≈ [2.0, 4.0] + Test.@test interp(0.5) ≈ [0.5, 1.0] + end + + Test.@testset "ctinterpolate_constant - basic piecewise constant" begin + x = [0.0, 1.0, 2.0] + f = [1.0, 2.0, 3.0] + interp = Interpolation.ctinterpolate_constant(x, f) + Test.@test interp(0.0) ≈ 1.0 + Test.@test interp(1.0) ≈ 2.0 + Test.@test interp(2.0) ≈ 3.0 + Test.@test interp(0.5) ≈ 1.0 + Test.@test interp(0.9) ≈ 1.0 + Test.@test interp(1.5) ≈ 2.0 + end + + Test.@testset "ctinterpolate_constant - extrapolation" begin + x = [0.0, 1.0, 2.0] + f = [1.0, 2.0, 3.0] + interp = Interpolation.ctinterpolate_constant(x, f) + Test.@test interp(-1.0) ≈ 1.0 + Test.@test interp(-0.5) ≈ 1.0 + Test.@test interp(3.0) ≈ 3.0 + Test.@test interp(5.0) ≈ 3.0 + end + + Test.@testset "ctinterpolate_constant - control-like data" begin + x = [0.0, 0.1, 0.5, 1.0] + f = [5.0, 4.0, 3.0, 2.0] + interp = Interpolation.ctinterpolate_constant(x, f) + Test.@test interp(0.0) ≈ 5.0 + Test.@test interp(0.05) ≈ 5.0 + Test.@test interp(0.1) ≈ 4.0 + Test.@test interp(0.3) ≈ 4.0 + Test.@test interp(0.5) ≈ 3.0 + Test.@test interp(0.75) ≈ 3.0 + Test.@test interp(1.0) ≈ 2.0 + end + + Test.@testset "ctinterpolate_constant - vector values" begin + x = [0.0, 1.0, 2.0] + f = [[1.0, 2.0], [3.0, 4.0], [5.0, 6.0]] + interp = Interpolation.ctinterpolate_constant(x, f) + Test.@test interp(0.0) ≈ [1.0, 2.0] + Test.@test interp(1.0) ≈ [3.0, 4.0] + Test.@test interp(2.0) ≈ [5.0, 6.0] + Test.@test interp(0.5) ≈ [1.0, 2.0] + Test.@test interp(1.5) ≈ [3.0, 4.0] + end + + Test.@testset "ctinterpolate_constant - constant function" begin + x = [0.0, 1.0, 2.0, 3.0] + f = [7.0, 7.0, 7.0, 7.0] + interp = Interpolation.ctinterpolate_constant(x, f) + Test.@test interp(0.5) ≈ 7.0 + Test.@test interp(1.5) ≈ 7.0 + Test.@test interp(2.5) ≈ 7.0 + Test.@test interp(-1.0) ≈ 7.0 + Test.@test interp(5.0) ≈ 7.0 + end + + end + return nothing +end + +end # module TestInterpolation + +test_interpolation() = TestInterpolation.test_interpolation() From d6d98bd3dad99c5c51218f479fadc791c61b65fb Mon Sep 17 00:00:00 2001 From: Olivier Cots Date: Thu, 11 Jun 2026 23:25:58 +0200 Subject: [PATCH 4/5] Refactor Interpolation to use typed Interpolant structs - Added Interpolant{Linear} and Interpolant{Constant} types - Implemented call methods for each type - Added display.jl with custom show methods - Updated docs to include types.jl and display.jl - Added type stability tests for interpolant calls --- docs/api_reference.jl | 4 +- src/Interpolation/Interpolation.jl | 5 +- src/Interpolation/ctinterpolate.jl | 79 ++++++++++--------- .../suite/interpolation/test_interpolation.jl | 20 +++++ 4 files changed, 70 insertions(+), 38 deletions(-) diff --git a/docs/api_reference.jl b/docs/api_reference.jl index d8a8d4a6..e9b850aa 100644 --- a/docs/api_reference.jl +++ b/docs/api_reference.jl @@ -39,12 +39,14 @@ function generate_api_reference(src_dir::String) primary_modules=[ CTBase.Interpolation => src( joinpath("Interpolation", "Interpolation.jl"), + joinpath("Interpolation", "types.jl"), joinpath("Interpolation", "ctinterpolate.jl"), + joinpath("Interpolation", "display.jl"), ), ], exclude=EXCLUDE_SYMBOLS, public=true, - private=false, + private=true, title="Interpolation", title_in_menu="Interpolation", filename="interpolation", diff --git a/src/Interpolation/Interpolation.jl b/src/Interpolation/Interpolation.jl index 03d0fc79..42607b3e 100644 --- a/src/Interpolation/Interpolation.jl +++ b/src/Interpolation/Interpolation.jl @@ -9,10 +9,13 @@ Interpolation utilities for the Control Toolbox (CT) ecosystem. """ module Interpolation -import DocStringExtensions: TYPEDSIGNATURES +import DocStringExtensions: TYPEDEF, TYPEDSIGNATURES +include(joinpath(@__DIR__, "types.jl")) include(joinpath(@__DIR__, "ctinterpolate.jl")) +include(joinpath(@__DIR__, "display.jl")) export ctinterpolate, ctinterpolate_constant +export Interpolant, LinearInterpolant, ConstantInterpolant end # module Interpolation diff --git a/src/Interpolation/ctinterpolate.jl b/src/Interpolation/ctinterpolate.jl index 35c4a749..e67e69b0 100644 --- a/src/Interpolation/ctinterpolate.jl +++ b/src/Interpolation/ctinterpolate.jl @@ -1,17 +1,52 @@ +# ============================================================================= +# Call methods +# ============================================================================= + +# Linear interpolation with flat extrapolation outside [x[1], x[end]]. +function (interp::Interpolant{Linear})(t) + x, f = interp.x, interp.f + if t < x[1] + return f[1] + elseif t >= x[end] + return f[end] + else + i = searchsortedlast(x, t) + i == length(x) && return f[end] + α = (t - x[i]) / (x[i + 1] - x[i]) + return f[i] + α * (f[i + 1] - f[i]) + end +end + +# Right-continuous piecewise-constant (steppost) interpolation. +function (interp::Interpolant{Constant})(t) + x, f = interp.x, interp.f + if t < x[1] + return f[1] + elseif t >= x[end] + return f[end] + else + return f[searchsortedlast(x, t)] + end +end + +# ============================================================================= +# Factories (public API) +# ============================================================================= + """ $(TYPEDSIGNATURES) Return a linear interpolation function for the data `f` defined at points `x`. -This function creates a one-dimensional linear interpolant with flat extrapolation -beyond the bounds of `x` (returns `f[1]` for `t < x[1]` and `f[end]` for `t >= x[end]`). +This creates a one-dimensional linear [`Interpolant`](@ref) with flat extrapolation beyond +the bounds of `x` (returns `f[1]` for `t < x[1]` and `f[end]` for `t >= x[end]`). # Arguments - `x`: A vector of points at which the values `f` are defined. - `f`: A vector of values to interpolate. # Returns -A callable interpolation object that can be evaluated at new points. +A callable [`LinearInterpolant`](@ref) that can be evaluated at new points. # Example ```julia-repl @@ -23,31 +58,15 @@ julia> interp(-1.0) # Returns 1.0 (flat extrapolation) julia> interp(3.0) # Returns 3.0 (flat extrapolation) ``` """ -function ctinterpolate(x, f) - function linear_interp(t) - if t < x[1] - return f[1] - elseif t >= x[end] - return f[end] - else - i = searchsortedlast(x, t) - if i == length(x) - return f[end] - end - α = (t - x[i]) / (x[i + 1] - x[i]) - return f[i] + α * (f[i + 1] - f[i]) - end - end - return linear_interp -end +ctinterpolate(x, f) = Interpolant{Linear}(x, f) """ $(TYPEDSIGNATURES) Return a piecewise-constant interpolation function for the data `f` defined at points `x`. -This function creates a right-continuous piecewise constant interpolant: -the value at knot `x[i]` is held constant on the interval `[x[i], x[i+1})`. +This creates a right-continuous piecewise-constant [`Interpolant`](@ref): the value at knot +`x[i]` is held constant on the interval `[x[i], x[i+1})`. This implements the standard steppost behavior for optimal control: - `u(t_i) = u_i` (value at the knot) @@ -59,7 +78,7 @@ This implements the standard steppost behavior for optimal control: - `f`: A vector of values to interpolate. # Returns -A callable interpolation object that can be evaluated at new points. +A callable [`ConstantInterpolant`](@ref) that can be evaluated at new points. # Example ```julia-repl @@ -72,16 +91,4 @@ julia> interp(1.0) # Returns 2.0 (value at x[2], right-continuous) julia> interp(1.5) # Returns 2.0 (held from x[2] on [1.0, 2.0)) ``` """ -function ctinterpolate_constant(x, f) - function steppost_interp(t) - if t < x[1] - return f[1] - elseif t >= x[end] - return f[end] - else - i = searchsortedlast(x, t) - return f[i] - end - end - return steppost_interp -end +ctinterpolate_constant(x, f) = Interpolant{Constant}(x, f) diff --git a/test/suite/interpolation/test_interpolation.jl b/test/suite/interpolation/test_interpolation.jl index 96544431..6b3438c8 100644 --- a/test/suite/interpolation/test_interpolation.jl +++ b/test/suite/interpolation/test_interpolation.jl @@ -113,6 +113,26 @@ function test_interpolation() Test.@test interp(5.0) ≈ 7.0 end + Test.@testset "Interpolant type & show" begin + interp = Interpolation.ctinterpolate([0.0, 1.0, 2.0], [0.0, 1.0, 0.0]) + Test.@test interp isa Function + Test.@test interp isa Interpolation.Interpolant + Test.@test interp isa Interpolation.LinearInterpolant + Test.@test Interpolation.method(interp) === Interpolation.Linear + Test.@test occursin("linear", sprint(show, interp)) + Test.@test occursin("3 nodes", sprint(show, interp)) + + c = Interpolation.ctinterpolate_constant([0.0, 1.0, 2.0], [1.0, 2.0, 3.0]) + Test.@test c isa Function + Test.@test c isa Interpolation.ConstantInterpolant + Test.@test Interpolation.method(c) === Interpolation.Constant + Test.@test occursin("constant", sprint(show, c)) + + # type stability of the call (philosophy: @inferred on the call) + Test.@test_nowarn Test.@inferred interp(0.5) + Test.@test_nowarn Test.@inferred c(0.5) + end + end return nothing end From d598d0fc52ad3ca055e8c37765d19587b87007d7 Mon Sep 17 00:00:00 2001 From: Olivier Cots Date: Thu, 11 Jun 2026 23:26:22 +0200 Subject: [PATCH 5/5] Add Interpolation types and display files --- src/Interpolation/display.jl | 24 +++++++++++++ src/Interpolation/types.jl | 68 ++++++++++++++++++++++++++++++++++++ 2 files changed, 92 insertions(+) create mode 100644 src/Interpolation/display.jl create mode 100644 src/Interpolation/types.jl diff --git a/src/Interpolation/display.jl b/src/Interpolation/display.jl new file mode 100644 index 00000000..42ad9a7c --- /dev/null +++ b/src/Interpolation/display.jl @@ -0,0 +1,24 @@ +""" +$(TYPEDSIGNATURES) + +Display a compact representation of an [`Interpolant`](@ref). + +# Example +```julia-repl +julia> ctinterpolate([0.0, 1.0, 2.0], [0.0, 1.0, 0.0]) +Interpolant (linear): 3 nodes +``` +""" +function Base.show(io::IO, interp::Interpolant{M}) where {M} + label = M === Linear ? "linear" : "piecewise-constant" + print(io, "Interpolant ($label): $(length(interp.x)) nodes") +end + +""" +$(TYPEDSIGNATURES) + +Display an [`Interpolant`](@ref) in the REPL with the same format as `Base.show(io, interp)`. +""" +function Base.show(io::IO, ::MIME"text/plain", interp::Interpolant) + return show(io, interp) +end diff --git a/src/Interpolation/types.jl b/src/Interpolation/types.jl new file mode 100644 index 00000000..5db421e8 --- /dev/null +++ b/src/Interpolation/types.jl @@ -0,0 +1,68 @@ +""" +$(TYPEDEF) + +Abstract interpolation method (trait). + +Concrete subtypes (`Linear`, `Constant`) are singleton trait values carried as a type +parameter of [`Interpolant`](@ref), so the call method is resolved statically. +""" +abstract type AbstractInterpolation end + +""" +$(TYPEDEF) + +Linear interpolation method (with flat extrapolation outside the node range). +""" +struct Linear <: AbstractInterpolation end + +""" +$(TYPEDEF) + +Piecewise-constant (right-continuous steppost) interpolation method. +""" +struct Constant <: AbstractInterpolation end + +""" +$(TYPEDEF) + +Callable interpolant of method `M` over nodes `x` with values `f`. + +`Interpolant` subtypes `Function`: an instance `interp` is evaluated as `interp(t)`. The +method `M` (a subtype of [`AbstractInterpolation`](@ref)) is a compile-time trait parameter +that selects the evaluation rule, so the call is type-stable. + +# Fields +- `x::TX`: nodes at which the values are defined. +- `f::TF`: values to interpolate. + +# Construction +Build instances through the factories [`ctinterpolate`](@ref) (linear) and +[`ctinterpolate_constant`](@ref) (piecewise-constant), or directly via the aliases +[`LinearInterpolant`](@ref) / [`ConstantInterpolant`](@ref). +""" +struct Interpolant{M<:AbstractInterpolation,TX,TF} <: Function + x::TX + f::TF +end + +# Outer constructor: infer TX, TF from the arguments. +function Interpolant{M}(x::TX, f::TF) where {M<:AbstractInterpolation,TX,TF} + return Interpolant{M,TX,TF}(x, f) +end + +""" +Alias for a linear [`Interpolant`](@ref). +""" +const LinearInterpolant = Interpolant{Linear} + +""" +Alias for a piecewise-constant [`Interpolant`](@ref). +""" +const ConstantInterpolant = Interpolant{Constant} + +""" +$(TYPEDSIGNATURES) + +Return the interpolation method (a subtype of [`AbstractInterpolation`](@ref)) of `interp`. +""" +method(::Interpolant{M}) where {M} = M