Search code examples
pythonscipyquaternionsodeintrigid-bodies

Unable to pass multiple arguments to the odeint solver


I am trying to numerically integrate the following equation where the omega(3x1) vector and theta are (4X1) vector and W(w) is 4x4 skew-symmetric matrix theta equation

def skew_omega(w):
    """skew symmetric omega (4x4)"""
        print(w)
        return np.array([[0, -w[0], -w[1], -w[2]],
                         [w[0], 0, w[2], -w[1]],
                         [w[1], -w[2], 0, w[0]],
                         [w[2], w[1], -w[0], 0]])


def theta_model(euler_param, w, t):
"""Ode for theta"""
    print(w)
    theta = 0.5*np.matmul(skew_omega(w),euler_param)
    
    return theta

#initial conditions
t = np.array([0, 1])
theta0 = np.array([1, 0, 0, 0])
omega0 = np.array([0.1, 10, 0])

#solve ode    
theta = odeint(theta_model, theta0, t, args = (omega0,))

0.0

0.0

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-72-1324ec16154f> in <module>
----> 1 theta = odeint(theta_model, theta0, t, args = (om,))

W:\Anaconda\lib\site-packages\scipy\integrate\odepack.py in odeint(func, y0, t, args, Dfun, col_deriv, full_output, ml, mu, rtol, atol, tcrit, h0, hmax, hmin, ixpr, mxstep, mxhnil, mxordn, mxords, printmessg, tfirst)
    243                              full_output, rtol, atol, tcrit, h0, hmax, hmin,
    244                              ixpr, mxstep, mxhnil, mxordn, mxords,
--> 245                              int(bool(tfirst)))
    246     if output[-1] < 0:
    247         warning_msg = _msgs[output[-1]] + " Run with full_output = 1 to get quantitative information."

<ipython-input-55-c581c2c300b8> in theta_model(euler_param, w, t)
      1 def theta_model(euler_param, w, t):
      2     print(w)
----> 3     theta = 0.5*np.matmul(skew_omega(w),euler_param)
      4 
      5     return theta

<ipython-input-49-bb06a7f0e258> in skew_omega(w)
      1 def skew_omega(w):
      2     print(w)
----> 3     return np.array([[0, -w[0], -w[1], -w[2]],
      4                      [w[0], 0, w[2], -w[1]],
      5                      [w[1], -w[2], 0, w[0]],

TypeError: 'float' object is not subscriptable

The printed zeros((the print statements in the functions) indicate that the odeint solver is not taking in the inital value- omega0. Instead it gives a 0.0 float as a place holder and hence the TypeError is thrown.

Why the ode solver not taking omega0 as input?


Solution

  • From odeint documentation

    Parameters
    ----------
    func : callable(y, t, ...) or callable(t, y, ...)
    Computes the derivative of y at t.

    Your model function

    def theta_model(euler_param, w, t):
    

    Looks like you have to exchange t and w in the argument list.