Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Documentation for Handling Trigonometric Functions in SINDy-PI Library Functions #498

Open
ibrahimw1 opened this issue May 8, 2024 · 1 comment

Comments

@ibrahimw1
Copy link

Firstly, thank you for making this package.
I'm currently working on a research project involving the application of SINDy-PI from the PySINDy library. My specific focus is on modeling a double pendulum system, which involves equations with trigonometric functions such as sine and cosine.

The current documentation for PySINDy, especially higher order ODEs, lacks clear guidance for handling library functions involving trigonometric operations like sine and cosine(also lacking general documentation for trig equations like double pendulum). While the provided examples and documentation cover high-order ODEs, they do not address scenarios where ODEs involve trigonometric terms(eg pendulum on cart and double pendulum).

I am trying to replicate the results of the double pendulum problem outlined in the SINDy-PI paper "SINDy-PI: A Robust Algorithm for Parallel Implicit Sparse Identification of Nonlinear Dynamics", where the equation that was found by SINDY-PI was:

image

So it would be very helpful if you could provide library functions for this problem. Please check out the code I have so far and provide any feedback:

import numpy as np
import math
from matplotlib import pyplot as plt
import pysindy as ps
from pysindy.differentiation import FiniteDifference
import sympy as sp
from scipy.integrate import solve_ivp
from pysindy.utils.odes import double_pendulum

dt = 0.001
x0_train = [np.pi + 1.2, np.pi - 0.6, 0, 0]
t_train = np.arange(0,10,dt)
t_train_span = (t_train[0], t_train[-1])
x_train = solve_ivp(double_pendulum, t_train_span, x0_train, t_eval=t_train).y.T

library_functions = []
library_function_names = []

# Create the PDELibrary
sindy_library = ps.PDELibrary(
    library_functions=library_functions,
    temporal_grid=t_train,
    function_names=library_function_names,
    include_bias=True,
    implicit_terms=True,
    derivative_order=2
)


lib = sindy_library.fit(x_train)
lib.transform(x_train)

print("With function names: ")
print(lib.get_feature_names(), "\n")


sindy_opt = ps.SINDyPI(
    threshold=0.0001,
    tol=1e-5,
    thresholder="l1",
    max_iter=6000,
)

model = ps.SINDy(
    optimizer=sindy_opt,
    feature_library=sindy_library,
)

model.fit(x_train, t=t_train)
model.print()

Additional info - when i include library functions such Polynomial tensored with Fourier, the model ends up outputting 100-200 equations, none of them giving me the desired output from the original SINDy-PI paper. If i don't include any library functions, most coefficients become nan(0).

@Jacob-Stevens-Haas
Copy link
Collaborator

Jacob-Stevens-Haas commented May 8, 2024

Hey, thanks for using the package! One of the challenges with SINDy-PI is that the number of models it searches grows combinatorially. The other issue is that there's often multiple ways of representing complex problems in pysindy. My guess is that the result in this paper required a lot of engineering.

That said, it's best if you can aggressively limit the number of features first to make sure you data generation is correct, then work more generally. Here's some suggestions:

  1. Add a control variable equal to $\phi_1 - \phi_2$
  2. Put them in a generalized library. Use inputs_per_library to make sure that they only act on the control variable.
  3. Try solving it without the polynomial library by explicitly calculating $\dot \phi_1$ and $\dot\phi_2$. You can then treat $\phi_1$ and $\phi_2$ as control variables.
  4. print lib.get_feature_names() to know what features to expect.

Admittedly, this is the part of the library that I understand the least.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants