Hello Blender Community,

I am new to the software and I am looking for some advice on the following problem.

Suppose that we are given the position and orientation of a rigid object in Euclidean space $(x,y,z)$ which are represented by the coordinates $(\boldsymbol{r},\Lambda)$ respectively. Further, let us assume that the evolution of the position and orientation is governed by an ordinary differential equation of the form $$(\dot{\boldsymbol{r}},\dot{\Lambda}) = f(\boldsymbol{r},\Lambda) $$ where dot's refer to time derivatives $\frac{d}{dt}$.

We then can solve this ODE numerically using, say, the Runge-Kutta algorithm. This then gives us a time indexed sequence of vectors and matrices, namely,
$$\boldsymbol{r}_i = \begin{bmatrix}x_i\\y_i \\ z_i \end{bmatrix} $$ for the object's position and

$$\Lambda_i = \begin{bmatrix}\Lambda_{11}^i& \Lambda_{12}^i & \Lambda_{13}^i\\\Lambda_{21}^i& \Lambda_{22}^i & \Lambda_{23}^i\\ \Lambda_{31}^i& \Lambda_{32}^i & \Lambda_{33}^i \end{bmatrix} $$ for it's orientation. If one were to take enough points $i$, then this sequence approximates the trajectory of the object (its position and orientation in space) by applying these sets of transformations to a chosen $(x,y,z)$ coordinate system.

I would like to animate the trajectory but I need some help/examples as to how to implement. I think a possible solution is as follows:

Step 1: Either solve the ODE in Python and import an array of solution points as a csv file OR set up the problem using C# (I think this is what Blender uses though I am at loss as to how to interface with it), implement Runga-Kutta and get an time indexed array of state-points.

Step 2: Interpolate the missing points to obtain a smooth trajectory.

Step 3: Fix a coordinate system to the animated object

Step 4: Somehow apply coordinate transformations $(\boldsymbol{r},\Lambda)$ to the object

Step 5: run the simulation from which we obtain an animation the trajectory of a rigid body in space.

I know how to solve ODE's in Python or MATLAB (Step 1) but I am at a loss as to how to complete steps 2-5. If anyone has any expertise with such a problem or can suggest any reading material/resources that would be greatly appreciated. Thank you.

Let me give an example. Let us consider the motion of a free rigid body. Euler tells us that the evolution of the body's motion is governed by $$\dot{\Omega} = \mathbb{I}\Omega\times\Omega$$ $$\dot{\Lambda} = \widehat{\Omega}\Lambda $$ where $\widehat{\Omega} = \Omega\times = \epsilon_{ijk}\Omega_k$, the antisymmetic tensor. To solve this equation we use the following routine in python:

from scipy.integrate import odeint
import numpy as np
# define the hat map
# takes in 3 vector and returns anti symmetric matrix
def hat(x):
    hat = np.array([[0,-x[2],x[1],[x[2],0,-x[0]],[-x[1],x[0],0]])
    return hat

# define vector field
# x[0] - x[2] Omega1-3
# x[3] - x[11] Lambda
# I moment of inertia

def f(x,t,I):
    dfdt = np.zeros(12)
    Om1 = x[0]
    Om2 = x[1]
    Om3 = x[2]
    La11 = x[3]
    La12 = x[4]
    La13 = x[5]
    La21 = x[6]
    La22 = x[7]
    La23 = x[8]
    La31 = x[9]
    La32 = x[10]
    La33 = x[11]
    Om = np.array([Om1,Om2,Om3])
    La = np.array([[La11,La12,La13],[La21,La22,La23],[La31,La32,La33]])
    Pi = np.matmult(I,Om)
    X = np.cross(Pi,Om)
    # Dynamics
    dfdt[0:4] = X
    # Reconstruction equations
    hatOm = hat(Om)
    Y = np.matmul(hatOm,La)
    Ry = np.reshape(Y,(9))
    dfdt[4:12] = Ry
    return dfdt

# simulate rigid body
# time interval
tspan = np.linspace(0,1,100)
# initial conditions
IC = np.zeros(12)
IC[0] = 1
# moment of inertia
I = np.array([[1,0,0],[0,1,0],[0,0,1]])
x = odeint(f,IC,tspan,args=(I))
Omega1 = x[0:,0]
Omega2 = x[0:,1]
Omega3 = x[0:,2]
La11 = x[0:,3]
La12 = x[0:,4]
La13 = x[0:,5]
La21 = x[0:,6]
La22 = x[0:,7]
La23 = x[0:,8]
La31 = x[0:,9]
La32 = x[0:,10]
La33 = x[0:,11]

My question is how to attach these simulated values to an object to be animated. I have tried the following code, yet the animation it produces does not seem realistic. I am not sure where the issue is.

# lets look at the motion in velocity space
import bpy
ob = bpy.data.objects["Sphere"]
start_position = (1,0,0)
frame_num = 0
t0 = len(tspan)
for i in range(t0-1):
    ob.location = (Omega1[i],Omega2[i],Omega3[i])
    ob.keyframe_insert(data_path='location',index = -1)
    frame_num +=20
  • $\begingroup$ Good question, the first part of your code has errors. Could you fix them please? Or even better, could your replace it with a simpler equation? $\endgroup$
    – Leander
    Commented Oct 28, 2019 at 16:07
  • $\begingroup$ What I don't understand is, that you have timesteps of 0.01 but jump 20 frames when settings keyframes. $\endgroup$
    – Leander
    Commented Oct 28, 2019 at 16:08
  • $\begingroup$ As for the motion path, simple create a new object and populate it with your samples, then connect the vertices. Bmesh. $\endgroup$
    – Leander
    Commented Oct 28, 2019 at 16:12


You must log in to answer this question.

Browse other questions tagged .