import numpy as np

from scipy.integrate import solve$\_$ivp

# Define the EP differential kinematic equations

def ep$\_$deriv$($t$,$ ep$,$ omega$)$:

q$0,$ q$1,$ q$2,$ q$3=$ ep

q$\_$dot $=$ np$.$zeros$(4)$

q$\_$dot$[0]=-0\mathrm{.}5*($q$1*$ omega$[0]+$ q$2*$ omega$[1]+$ q$3*$ omega$[2])$

q$\_$dot$[1]=0\mathrm{.}5*($q$0*$ omega$[0]-$ q$3*$ omega$[1]+$ q$2*$ omega$[2])$

q$\_$dot$[2]=0\mathrm{.}5*($q$3*$ omega$[0]+$ q$0*$ omega$[1]-$ q$1*$ omega$[2])$

q$\_$dot$[3]=0\mathrm{.}5*(-$q$2*$ omega$[0]+$ q$1*$ omega$[1]+$ q$0*$ omega$[2])$

return q$\_$dot

# Define the body angular velocity as a function of time

def angular$\_$velocity$($t$)$:

return np$.$array$([$np$.$sin$(0\mathrm{.}1*$ t$),0\mathrm{.}01,$ np$.$cos$(0\mathrm{.}1*$ t$)])$

# Initial EP vector

initial$\_$ep $=$ np$.$array$([0\mathrm{.}408248,0.,0\mathrm{.}408248,0\mathrm{.}816497])$

# Time span for integration

t$\_$span $=(0,42)$

# Integrate the EP differential kinematic equations

solution $=$ solve$\_$ivp$($

fun$=$lambda t$,$ ep: ep$\_$deriv$($t$,$ ep$,$ angular$\_$velocity$($t$)),$

t$\_$span$=$t$\_$span,

y$0=$initial$\_$ep$,$

method$=\text{'}$RK$45\text{'},$

dense$\_$output$=$True,

$)$

# Evaluate the solution at t$=42$ seconds

result$\_$at$\_42\_$seconds $=$ solution.sol$(42)$

# Calculate the norm of the EP vector component at $42$ seconds

norm$\_$ep $=$ np$.$linalg.norm$($result$\_$at$\_42\_$seconds$[1$:$])$

print$($f$"$The norm of the EP vector component at $42$ seconds is: $\{$norm$\_$ep$\}")$