JuliaGNI Integration

The package is built on the JuliaGNI stack and composes with it in both directions: an integrated solution can be turned into training data, and an identified system can be turned back into a problem an integrator accepts.

The shape of the API

Identification mirrors integration. Where GeometricIntegrators has

solution = integrate(problem, method)

this package has

result = identify(problem, method)

TrainingData and TrajectoryData are GeometricBase.AbstractProblems, and the methods are GeometricBase.AbstractMethods, so the ecosystem's accessors work on them: nsamples, statedimension, timestep, datatype, arrtype, parameters, functions.

SparseIdentification.identify — Function

identify(problem::IdentificationProblem, method::SparsificationMethod; kwargs...)

Identify the governing equations of problem with method.

This mirrors GeometricIntegrators' integrate(problem, method): the problem carries the data, the method carries the basis and the sparsification parameters, and the result carries the identified coefficients.

Which methods accept which problems is not interchangeable, because the two formulations need different data:

problem	method	fits
`TrainingData`	`SINDy`	a vector field, against measured derivatives
`TrajectoryData`	`HamiltonianSINDy`	a Hamiltonian, against the flow map

Applying a method to a problem it does not accept raises an ArgumentError naming both, rather than a MethodError from several layers down.

Examples

julia> A = [-0.1 2.0; -2.0 -0.1];

julia> x = randn(2, 400);

julia> result = identify(TrainingData(x, A * x), SINDy(CompoundBasis(polyorder = 3); λ = 0.05));

julia> nterms(result)
4

source

SparseIdentification.IdentificationProblem — Type

IdentificationProblem

Supertype of the data an identification method is applied to.

An identification problem plays the role that a GeometricEquations problem plays for an integrator: it is the thing a method is applied to. identify is the verb, so that

identify(problem, method)

reads as GeometricIntegrators' integrate(problem, method) does.

Two concrete problems, distinguished by what the data means rather than by its shape: TrainingData for matching a vector field, TrajectoryData for matching a flow map.

source

SparseIdentification.SparsificationMethod — Type

SparsificationMethod

Supertype of the identification methods, e.g. SINDy and HamiltonianSINDy.

A method carries the basis it searches and the parameters of its sparsification. It is applied to an IdentificationProblem with identify, mirroring the way a GeometricIntegrators method is applied to a problem with integrate.

Traits

SparsificationMethod extends GeometricBase.AbstractMethod, and the ecosystem's trait functions are answered where they mean something for a regression method:

issymplectic and isenergypreserving report whether the identified model is symplectic and energy-preserving by construction. This is the substantive question for a structure-preserving method, and the reason HamiltonianSINDy exists.
isexplicit and isimplicit report whether the regression is a direct solve or goes through an implicit integrator.
name, description and reference identify the method and its source.

issymmetric, isstifflyaccurate and order describe a Runge–Kutta tableau and have no meaning for a regression, so they are left as missing. GeometricBase.isAbstractMethod therefore returns false for these methods, which is correct: it is a conformance check for integrator methods.

source

Closing the loop

The full round trip — integrate a known system, identify it from the solution, and integrate the identified system:

using SparseIdentification
using GeometricIntegrators
using GeometricProblems.HarmonicOscillator

# 1. a reference problem, integrated
reference = hodeproblem()
solution  = integrate(reference, ImplicitMidpoint())

# 2. its solution as training data — for a Hamiltonian system the state is z = (q, p)
data = TrajectoryData(solution)
(nsamples(data), statedimension(data))

(10, 2)

# 3. identify
method = HamiltonianSINDy(λ = 0.01, integrator_timestep = timestep(data), polyorder = 2)
result = identify(data, method)

Hamiltonian SINDy result: 2 of 5 coefficients retained, 1 degrees of freedom

# 4. back to a problem, and integrate it
identified = HODEProblem(result, (0.0, 1.0), 0.01, [0.5], [0.0])
idsolution = integrate(identified, ImplicitMidpoint())

typeof(idsolution).name.name

:GeometricSolution

Because the identified problem carries the identified Hamiltonian as an invariant, the energy behaviour of the integrated solution can be checked with GeometricSolutions' own diagnostics rather than by hand.

SparseIdentification.TrainingData — Type

TrainingData(x, ẋ)
TrainingData(solution::GeometricSolution)

States and the corresponding time derivatives, for methods that match a vector field.

x and ẋ are either matrices whose columns are snapshots, or vectors of state vectors. Both must describe the same number of snapshots of the same dimension; that is checked at construction, because the alternative is a shape error surfacing much later as a wrong answer.

Examples

julia> data = TrainingData(randn(2, 20), randn(2, 20));

julia> nsamples(data), statedimension(data)
(20, 2)

source

SparseIdentification.TrajectoryData — Type

TrajectoryData(x, y, Δt)
TrajectoryData(solution::GeometricSolution)

Consecutive states, for methods that match a flow map: y[j] is the state one step of size Δt after x[j].

Use this where the derivatives are not available and only sampled trajectories are — matching the flow map avoids differentiating the data, at the cost of a nonlinear regression.

Examples

julia> data = TrajectoryData([randn(2) for _ in 1:8], [randn(2) for _ in 1:8], 0.01);

julia> nsamples(data), statedimension(data), timestep(data)
(8, 2, 0.01)

source

Method traits

The ecosystem's trait functions are answered where they mean something for a regression method. The substantive ones are whether the identified model is structure-preserving:

using SparseIdentification: issymplectic, isenergypreserving

sindy = SINDy(CompoundBasis(polyorder = 2); λ = 0.05)
ham   = HamiltonianSINDy(λ = 0.05, polyorder = 2)

(issymplectic(sindy), issymplectic(ham))

(false, true)

That is the difference between the two methods stated as a property rather than as prose.

These six traits are not exported

isexplicit, isimplicit, issymmetric, issymplectic, isenergypreserving and isstifflyaccurate are extended here but reached qualified, as above.

GeometricIntegratorsBase defines its own generic functions of those six names rather than extending GeometricBase's stubs, and exports them. A session with both using SparseIdentification and using GeometricIntegrators would therefore see two different bindings under one name and resolve it to neither. Not exporting them is the same choice SimpleSolvers makes for status and isconverged, for the same reason.

issymmetric, isstifflyaccurate and order describe a Runge–Kutta tableau and have no meaning for a regression, so they are left as missing. GeometricBase.isAbstractMethod consequently returns false for these methods — which is correct, since it is a conformance check for integrator methods, not a general one.

What is deliberately not reused

GeometricBase.GeometricData wraps a NamedTuple of series tagged by system and data type. It is not used here: it carries no dimensions or sample count, has no outer constructor, and adding one from this package would be type piracy. TrainingData and TrajectoryData validate their shapes at construction and answer the ecosystem's accessors, which is what the wrapper would have been for.