Equations of Motion v1#
This notebook briefly states the equations of motion used in RocketPy, using symbolic algebra to prepare them for numerical integration. The result is used inside the Flight class to simulate a rocket flight.
Additional Requirements#
The code in this notebook requires sympy, Cython, Aesara and g++, which are not generally prerequisites to run RocketPy. You can install them with the next cells.
[ ]:
%pip install --upgrade sympy
[ ]:
%pip install --upgrade cython
[ ]:
%pip install --upgrade aesara
[ ]:
%conda install m2w64-toolchain
Importing libraries#
[4]:
import numpy as np
from sympy import *
from sympy.utilities.autowrap import autowrap, ufuncify, binary_function
from sympy.printing.aesaracode import aesara_function
from sympy.physics.vector import (
ReferenceFrame,
dynamicsymbols,
express,
time_derivative,
dot,
init_vprinting,
)
from sympy.physics.mechanics import inertia, msubs
# init_vprinting()
Example with Simple Harmonic Motion#
[5]:
x_t, y_t = dynamicsymbols("x y")
m1, m2, rho, k, t = symbols("m1 m2 rho k t")
x, xd, y, yd = symbols("x xd y yd")
xydd = solve(
[
Eq(
m1 * Derivative(x_t, t, 2)
+ rho * (Derivative(x_t, t, 1) - Derivative(y_t, t, 1))
+ k * (x_t - y_t),
0,
),
Eq(
m2 * Derivative(y_t, t, 2)
- rho * (Derivative(x_t, t, 1) - Derivative(y_t, t, 1))
- k * (x_t - y_t),
0,
),
],
[
Derivative(x_t, t, 2),
Derivative(y_t, t, 2),
],
)
print("x'':", xydd[Derivative(x_t, t, 2)])
print("y'':", xydd[Derivative(y_t, t, 2)])
xdd_clean = xydd[Derivative(x_t, t, 2)].subs(
{x_t: x, Derivative(x_t, t, 1): xd, y_t: y, Derivative(y_t, t, 1): yd}
)
ydd_clean = xydd[Derivative(y_t, t, 2)].subs(
{x_t: x, Derivative(x_t, t, 1): xd, y_t: y, Derivative(y_t, t, 1): yd}
)
print("x'':", xdd_clean)
print("y'':", ydd_clean)
xdd_function_eval_1 = lambdify((m1, m2, k, rho, x, xd, y, yd), xdd_clean)
xdd_function_eval_2 = ufuncify(
(m1, m2, k, rho, x, xd, y, yd), xdd_clean, backend="cython"
)
xdd_function_eval_3 = autowrap(xdd_clean, backend="cython")
xdd_function_eval_4 = aesara_function([m1, k, rho, x, xd, y, yd], [xdd_clean])
print(xdd_function_eval_1)
print(xdd_function_eval_2)
print(xdd_function_eval_3)
print(xdd_function_eval_4)
x'': -k*x(t)/m1 + k*y(t)/m1 - rho*Derivative(x(t), t)/m1 + rho*Derivative(y(t), t)/m1
y'': k*x(t)/m2 - k*y(t)/m2 + rho*Derivative(x(t), t)/m2 - rho*Derivative(y(t), t)/m2
x'': -k*x/m1 + k*y/m1 - rho*xd/m1 + rho*yd/m1
y'': k*x/m2 - k*y/m2 + rho*xd/m2 - rho*yd/m2
<function _lambdifygenerated at 0x000001C37945C040>
<built-in function autofunc_c>
<built-in function autofunc_c>
<aesara.compile.function.types.Function object at 0x000001C37939CA30>
[6]:
%timeit xdd_function_eval_1(1, 1, 1, 1, 1, 1, 1, 1)
args = (
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
np.array([1], dtype=np.double),
)
%timeit xdd_function_eval_2(*args)
%timeit xdd_function_eval_3(1, 1, 1, 1, 1, 1, 1)
%timeit xdd_function_eval_4(1, 1, 1, 1, 1, 1, 1)
711 ns ± 64.3 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
2.88 µs ± 125 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)
310 ns ± 7.97 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)
272 µs ± 41.8 µs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)
Six-Degree of Freedom Flight Phase#
[37]:
# Define state variables
x_t, y_t, z_t, q0_t, q1_t, q2_t, q3_t, w1_t, w2_t, w3_t = dynamicsymbols(
"x y z q_0 q_1 q_2 q_3 omega_1 omega_2 omega_3"
)
# Define state variable derivatives
vx_t, vy_t, vz_t = dynamicsymbols("x y z", 1)
ax_t, ay_t, az_t = dynamicsymbols("x y z", 3)
q0d_t, q1d_t, q2d_t, q3d_t = dynamicsymbols("q0 q1 q2 q3", 1)
w1d_t, w2d_t, w3d_t = dynamicsymbols("w1 w2 w3", 1)
# Define Reference Frames
## Inertial - Fixed to Earth
A = ReferenceFrame("A")
## Fixed to the rocket body (dry solid mass)
B = ReferenceFrame("B")
B.orient_quaternion(A, (q0_t, q1_t, q2_t, q3_t))
# Define parameters
m_t = dynamicsymbols("m", positive=True)
md_t = dynamicsymbols("m", 1)
v_t = vx_t * A.x + vy_t * A.y + vz_t * A.z
omega_t = w1_t * B.x + w2_t * B.y + w3_t * B.z
r_cm_scalar_t = dynamicsymbols("r_cm")
r_cm_t = r_cm_scalar_t * B.z
r_noz_scalar = symbols("r_noz")
r_noz = r_noz_scalar * B.z
A3_t = dynamicsymbols("A_3")
A_vec_t = -A3_t * B.z
Nx_t, Ny_t = dynamicsymbols("N_x, N_y")
N_t = Nx_t * B.x + Ny_t * B.y
T_scalar_t = dynamicsymbols("T")
T_t = T_scalar_t * B.z
g = symbols("g")
# Define translational equations of motion
trans_eq_lhs = (
m_t
* (
time_derivative(v_t, A)
+ time_derivative(omega_t, A).cross(r_cm_t)
+ omega_t.cross(
omega_t.cross(r_cm_t)
+ time_derivative(r_cm_t, B, 2)
+ 2 * omega_t.cross(time_derivative(r_cm_t, B))
)
)
).to_matrix(B)
trans_eq_rhs = (
T_t
- 2 * time_derivative(m_t, A) * time_derivative(r_cm_t, B)
+ 2 * omega_t.cross(time_derivative(m_t, A) * (r_noz - r_cm_t))
+ time_derivative(m_t, A, 2) * (r_noz - r_cm_t)
+ A_vec_t
+ N_t
- m_t * g * A.z
).to_matrix(B)
trans_eqs = [
Eq(trans_eq_lhs[0], trans_eq_rhs[0]),
Eq(trans_eq_lhs[1], trans_eq_rhs[1]),
Eq(trans_eq_lhs[2], trans_eq_rhs[2]),
]
# Parameters for rotational equations of motion
## Rocket Moment of Inertia
Ixx_t, Iyy_t, Izz_t, Ixy_t, Iyz_t, Izx_t = dynamicsymbols(
"I_xx I_yy I_zz I_xy I_yz I_zx"
)
I_t = inertia(B, Ixx_t, Iyy_t, Izz_t, Ixy_t, Iyz_t, Izx_t)
## Nozzle Gyration Tensor
S_noz_xx, S_noz_yy, S_noz_zz, S_noz_xy, S_noz_yz, S_noz_zx = symbols(
"S_noz_xx S_noz_yy S_noz_zz S_noz_xy S_noz_yz S_noz_zx"
)
S_noz = inertia(B, S_noz_xx, S_noz_yy, S_noz_zz, S_noz_xy, S_noz_yz, S_noz_zx)
## Moments
Mx_t, My_t, Mz_t = dynamicsymbols("M_x M_y M_z")
M_t = Mx_t * B.x + My_t * B.y + Mz_t * B.z
# Define rotational equations of motion
rot_eq_lhs = (
I_t.dot(time_derivative(omega_t, B))
+ omega_t.cross(I_t.dot(omega_t))
+ time_derivative(I_t, B).dot(omega_t)
+ m_t * r_cm_t.cross(time_derivative(v_t, A))
).to_matrix(B)
rot_eq_rhs = (
time_derivative(m_t, B) * S_noz.dot(omega_t) - r_cm_t.cross(m_t * g * (A.z)) + M_t
).to_matrix(B)
rot_eqs = [
Eq(rot_eq_lhs[0], rot_eq_rhs[0]),
Eq(rot_eq_lhs[1], rot_eq_rhs[1]),
Eq(rot_eq_lhs[2], rot_eq_rhs[2]),
]
eoms = trans_eqs + rot_eqs
# Clean dynamic variables
x, y, z, vx, vy, vz, q0, q1, q2, q3, w1, w2, w3 = symbols(
"x y z v_x v_y v_z q_0 q_1 q_2 q_3 omega_1 omega_2 omega_3"
)
ax, ay, az, q0d, q1d, q2d, q3d, w1d, w2d, w3d = symbols(
"a_x a_y a_z q_0d q_1d q_2d q_3d omega_1d omega_2d omega_3d"
)
m, md, mdd = symbols("m m_d m_dd")
r_cm_scalar, r_cmd, r_cmdd, r_noz_scalar, A3, Nx, Ny, T_scalar = symbols(
"r_cm r_cmd r_cmdd r_noz A_3 N_x N_y T"
)
Ixx, Iyy, Izz, Ixy, Iyz, Izx = symbols("I_xx I_yy I_zz I_xy I_yz I_zx")
Ixxd, Iyyd, Izzd, Ixyd, Iyzd, Izxd = symbols("I_xxd I_yyd I_zzd I_xyd I_yzd I_zxd")
Mx, My, Mz = symbols("M_x M_y M_z")
subs_dir = {
x_t: x,
y_t: y,
z_t: z,
vx_t: vx,
vy_t: vy,
vz_t: vz,
q0_t: q0,
q1_t: q1,
q2_t: q2,
q3_t: q3,
w1_t: w1,
w2_t: w2,
w3_t: w3,
Derivative(vx_t, t): ax,
Derivative(vy_t, t): ay,
Derivative(vz_t, t): az,
ax_t: ax,
ay_t: ay,
az_t: az,
q0d_t: q0d,
q1d_t: q1d,
q2d_t: q2d,
q3d_t: q3d,
Derivative(q0_t, t): q0d,
Derivative(q1_t, t): q1d,
Derivative(q2_t, t): q2d,
Derivative(q3_t, t): q3d,
w1d_t: w1d,
w2d_t: w2d,
w3d_t: w3d,
Derivative(w1_t, t): w1d,
Derivative(w2_t, t): w2d,
Derivative(w3_t, t): w3d,
m_t: m,
md_t: md,
Derivative(m_t, t): md,
Derivative(m_t, t, 2): mdd,
r_cm_scalar_t: r_cm_scalar,
r_noz_scalar: r_noz_scalar,
A3_t: A3,
Derivative(r_cm_scalar_t, t): r_cmd,
Derivative(r_cm_scalar_t, t, 2): r_cmdd,
Nx_t: Nx,
Ny_t: Ny,
T_scalar_t: T_scalar,
Ixx_t: Ixx,
Iyy_t: Iyy,
Izz_t: Izz,
Ixy_t: Ixy,
Iyz_t: Iyz,
Izx_t: Izx,
Derivative(Ixx_t, t): Ixxd,
Derivative(Iyy_t, t): Iyyd,
Derivative(Izz_t, t): Izzd,
Derivative(Ixy_t, t): Ixyd,
Derivative(Iyz_t, t): Iyzd,
Derivative(Izx_t, t): Izxd,
Mx_t: Mx,
My_t: My,
Mz_t: Mz,
}
eoms_cleaned = [msubs(eom, subs_dir) for eom in eoms]
# Convert EOMs to matrix form
eoms_cleaned_matrix_form = linear_eq_to_matrix(
eoms_cleaned, *[ax, ay, az, w1d, w2d, w3d]
)
AA = eoms_cleaned_matrix_form[0]
bb = eoms_cleaned_matrix_form[1]
# Compile AA and bb
AA_func_cython = autowrap(
AA,
backend="cython",
args=[m, r_cm_scalar, q0, q1, q2, q3, Ixx, Iyy, Izz, Ixy, Izx, Iyz],
tempdir="./compiled_funcs",
)
AA_func_aesara = aesara_function(
[m, r_cm_scalar, q0, q1, q2, q3, Ixx, Iyy, Izz, Ixy, Izx, Iyz],
[AA],
on_unused_input="warn",
)
bb_func_cython = autowrap(
bb,
backend="cython",
args=[
q0,
q1,
q2,
q3,
w1,
w2,
w3,
q0d,
q1d,
q2d,
q3d,
m,
md,
mdd,
r_cm_scalar,
r_cmd,
r_cmdd,
r_noz_scalar,
A3,
Nx,
Ny,
T_scalar,
Ixx,
Iyy,
Izz,
Ixy,
Iyz,
Izx,
Ixxd,
Iyyd,
Izzd,
Ixyd,
Iyzd,
Izxd,
Mx,
My,
Mz,
g,
S_noz_xx,
S_noz_yy,
S_noz_zz,
S_noz_xy,
S_noz_yz,
S_noz_zx,
],
tempdir="./compiled_funcs",
)
bb_func_aesara = aesara_function(
[
q0,
q1,
q2,
q3,
w1,
w2,
w3,
q0d,
q1d,
q2d,
q3d,
m,
md,
mdd,
r_cm_scalar,
r_cmd,
r_cmdd,
r_noz_scalar,
A3,
Nx,
Ny,
T_scalar,
Ixx,
Iyy,
Izz,
Ixy,
Iyz,
Izx,
Ixxd,
Iyyd,
Izzd,
Ixyd,
Iyzd,
Izxd,
Mx,
My,
Mz,
g,
S_noz_xx,
S_noz_yy,
S_noz_zz,
S_noz_xy,
S_noz_yz,
S_noz_zx,
],
[bb],
on_unused_input="warn",
)
# Solve for accelerations
# accelerations = linsolve((AA, bb))
# print(accelerations)
[38]:
AA
[38]:
$\displaystyle \left[\begin{matrix}m \left(q_{0}^{2} + q_{1}^{2} - q_{2}^{2} - q_{3}^{2}\right) & m \left(2 q_{0} q_{3} + 2 q_{1} q_{2}\right) & m \left(- 2 q_{0} q_{2} + 2 q_{1} q_{3}\right) & 0 & m r_{cm} & 0\\m \left(- 2 q_{0} q_{3} + 2 q_{1} q_{2}\right) & m \left(q_{0}^{2} - q_{1}^{2} + q_{2}^{2} - q_{3}^{2}\right) & m \left(2 q_{0} q_{1} + 2 q_{2} q_{3}\right) & - m r_{cm} & 0 & 0\\m \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) & m \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) & m \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) & 0 & 0 & 0\\- m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(- 2 q_{0} q_{2} + 2 q_{1} q_{3}\right) + m r_{cm} \left(2 q_{0} q_{3} + 2 q_{1} q_{2}\right) \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) & m r_{cm} \left(- 2 q_{0} q_{2} + 2 q_{1} q_{3}\right) \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) - m r_{cm} \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) \left(q_{0}^{2} + q_{1}^{2} - q_{2}^{2} - q_{3}^{2}\right) & m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(q_{0}^{2} + q_{1}^{2} - q_{2}^{2} - q_{3}^{2}\right) - m r_{cm} \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) \left(2 q_{0} q_{3} + 2 q_{1} q_{2}\right) & I_{xx} & I_{xy} & I_{zx}\\- m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(2 q_{0} q_{1} + 2 q_{2} q_{3}\right) + m r_{cm} \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) \left(q_{0}^{2} - q_{1}^{2} + q_{2}^{2} - q_{3}^{2}\right) & m r_{cm} \left(2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) - m r_{cm} \left(- 2 q_{0} q_{3} + 2 q_{1} q_{2}\right) \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) & m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(- 2 q_{0} q_{3} + 2 q_{1} q_{2}\right) - m r_{cm} \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) \left(q_{0}^{2} - q_{1}^{2} + q_{2}^{2} - q_{3}^{2}\right) & I_{xy} & I_{yy} & I_{yz}\\0 & 0 & 0 & I_{zx} & I_{yz} & I_{zz}\end{matrix}\right]$
[39]:
bb
[39]:
$\displaystyle \left[\begin{matrix}N_{x} - g m \left(- 2 q_{0} q_{2} + 2 q_{1} q_{3}\right) - m \left(\omega_{2} r_{cmdd} - \omega_{3} \left(- \omega_{1} r_{cm} - 2 \omega_{1} r_{cmd}\right) + r_{cm} \left(\omega_{1} \cdot \left(2 q_{0} q_{3d} - 2 q_{0d} q_{3} - 2 q_{1} q_{2d} + 2 q_{1d} q_{2}\right) - \omega_{3} \cdot \left(2 q_{0} q_{1d} - 2 q_{0d} q_{1} - 2 q_{2} q_{3d} + 2 q_{2d} q_{3}\right)\right)\right) + 2 m_{d} \omega_{2} \left(- r_{cm} + r_{noz}\right)\\N_{y} - g m \left(2 q_{0} q_{1} + 2 q_{2} q_{3}\right) - m \left(- \omega_{1} r_{cmdd} + \omega_{3} \left(\omega_{2} r_{cm} + 2 \omega_{2} r_{cmd}\right) - r_{cm} \left(- \omega_{2} \cdot \left(2 q_{0} q_{3d} - 2 q_{0d} q_{3} - 2 q_{1} q_{2d} + 2 q_{1d} q_{2}\right) + \omega_{3} \cdot \left(2 q_{0} q_{2d} - 2 q_{0d} q_{2} + 2 q_{1} q_{3d} - 2 q_{1d} q_{3}\right)\right)\right) - 2 m_{d} \omega_{1} \left(- r_{cm} + r_{noz}\right)\\- A_{3} + T - g m \left(q_{0}^{2} - q_{1}^{2} - q_{2}^{2} + q_{3}^{2}\right) - m \left(\omega_{1} \left(- \omega_{1} r_{cm} - 2 \omega_{1} r_{cmd}\right) - \omega_{2} \left(\omega_{2} r_{cm} + 2 \omega_{2} r_{cmd}\right)\right) - 2 m_{d} r_{cmd} + m_{dd} \left(- r_{cm} + r_{noz}\right)\\- I_{xxd} \omega_{1} - I_{xyd} \omega_{2} - I_{zxd} \omega_{3} + M_{x} - g m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(q_{0}^{2} + q_{1}^{2} - q_{2}^{2} - q_{3}^{2}\right) + g m r_{cm} \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) \left(2 q_{0} q_{3} + 2 q_{1} q_{2}\right) + m_{d} \left(S_{noz xx} \omega_{1} + S_{noz xy} \omega_{2} + S_{noz zx} \omega_{3}\right) - \omega_{2} \left(I_{yz} \omega_{2} + I_{zx} \omega_{1} + I_{zz} \omega_{3}\right) + \omega_{3} \left(I_{xy} \omega_{1} + I_{yy} \omega_{2} + I_{yz} \omega_{3}\right)\\- I_{xyd} \omega_{1} - I_{yyd} \omega_{2} - I_{yzd} \omega_{3} + M_{y} - g m r_{cm} \left(- 2 q_{0} q_{1} + 2 q_{2} q_{3}\right) \left(- 2 q_{0} q_{3} + 2 q_{1} q_{2}\right) + g m r_{cm} \left(2 q_{0} q_{2} + 2 q_{1} q_{3}\right) \left(q_{0}^{2} - q_{1}^{2} + q_{2}^{2} - q_{3}^{2}\right) + m_{d} \left(S_{noz xy} \omega_{1} + S_{noz yy} \omega_{2} + S_{noz yz} \omega_{3}\right) + \omega_{1} \left(I_{yz} \omega_{2} + I_{zx} \omega_{1} + I_{zz} \omega_{3}\right) - \omega_{3} \left(I_{xx} \omega_{1} + I_{xy} \omega_{2} + I_{zx} \omega_{3}\right)\\- I_{yzd} \omega_{2} - I_{zxd} \omega_{1} - I_{zzd} \omega_{3} + M_{z} + m_{d} \left(S_{noz yz} \omega_{2} + S_{noz zx} \omega_{1} + S_{noz zz} \omega_{3}\right) - \omega_{1} \left(I_{xy} \omega_{1} + I_{yy} \omega_{2} + I_{yz} \omega_{3}\right) + \omega_{2} \left(I_{xx} \omega_{1} + I_{xy} \omega_{2} + I_{zx} \omega_{3}\right)\end{matrix}\right]$
Testing Compiled Functions#
[35]:
import sys
sys.path.append("./compiled_funcs")
import wrapper_module_14 as AA
import wrapper_module_15 as bb
m = 20
r_cm_scalar = 0.1
q0 = 1
q1 = 0
q2 = 0
q3 = 0
Ixx = 1
Iyy = 1
Izz = 0.1
Ixy = 0
Izx = 0
Iyz = 0
w1 = 0
w2 = 0
w3 = 0
q0d = 0
q1d = 0
q2d = 0
q3d = 0
md = -1
mdd = 0
r_cmd = -0.01
r_cmdd = 0
r_noz_scalar = 1
A3 = 10
Nx = 0
Ny = 0
T_scalar = 1000
Ixxd = -0.1
Iyyd = -0.1
Izzd = -0.1
Ixyd = 0
Iyzd = 0
Izxd = 0
Mx = 0
My = 0
Mz = 0
g = -9.8
S_noz_xx = 1
S_noz_yy = 1
S_noz_zz = 1
S_noz_xy = 0
S_noz_yz = 0
S_noz_zx = 0
print(AA.autofunc_c(m, r_cm_scalar, q0, q1, q2, q3, Ixx, Iyy, Izz, Ixy, Izx, Iyz))
print(
bb.autofunc_c(
q0,
q1,
q2,
q3,
w1,
w2,
w3,
q0d,
q1d,
q2d,
q3d,
m,
md,
mdd,
r_cm_scalar,
r_cmd,
r_cmdd,
r_noz_scalar,
A3,
Nx,
Ny,
T_scalar,
Ixx,
Iyy,
Izz,
Ixy,
Iyz,
Izx,
Ixxd,
Iyyd,
Izzd,
Ixyd,
Iyzd,
Izxd,
Mx,
My,
Mz,
g,
S_noz_xx,
S_noz_yy,
S_noz_zz,
S_noz_xy,
S_noz_yz,
S_noz_zx,
)
)
[[20. 0. 0. 0. 2. 0. ]
[ 0. 20. 0. -2. 0. 0. ]
[ 0. 0. 20. 0. 0. 0. ]
[ 0. -2. 0. 1. 0. 0. ]
[ 2. 0. 0. 0. 1. 0. ]
[ 0. 0. 0. 0. 0. 0.1]]
[[ 0. ]
[ 0. ]
[1185.98]
[ 0. ]
[ 0. ]
[ 0. ]]