Motion in 3D.sagews

# Motion in 3D
# Author: Thomas Krainer (Penn State Altoona)
# Version: 02/11/2016
# Prior Version: 08/26/2011
#
# Changes: Norms changed from method .norm() to calculations by hand to fix
# issues with running the program.
#
# Based on "3-D motion and vector calculus" by Robert A. Beezer
# See http://wiki.sagemath.org/interact/calculus
#

t=var('t')
u,v=var('u,v')

output='<h1>Motion in 3-dimensional space</h1>'

output+='$$ \\begin{aligned} &\\textrm{Position:} &&\\vec{p}(t) = \\langle x(t),y(t),z(t) \\rangle'
output+=' \\quad \\textrm{for } t \in [t_{\\min},t_{\\max}] \\\ '
output+='&\\textrm{Velocity:} &&\\vec{v} = \\frac{d\\vec{p}}{dt} \\\ '
output+='&\\textrm{Speed:} &&v = |\\vec{v}| \\\ '
output+='&\\textrm{Acceleration:} &&\\vec{a} = \\frac{d\\vec{v}}{dt} = \\frac{d^2\\vec{p}}{dt^2} = a_{T}\\vec{T} + a_N\\vec{N} = \\left(\\frac{dv}{dt}\\right) \\vec{T} + \\left(\\kappa v^2\\right) \\vec{N} = \\left(\\frac{\\vec{v}\\cdot\\vec{a}}{v}\\right)\\vec{T} + \\left(\\frac{|\\vec{v}\\times\\vec{a}|}{v}\\right)\\vec{N} \\\ '
output+='&\\textrm{Unit tangent:} &&\\vec{T} = \\frac{\\vec{v}}{v} \\\ '
output+='&\\textrm{Unit normal:} &&\\vec{N} = \\frac{\\frac{d\\vec{T}}{dt}}{\\left|\\frac{d\\vec{T}}{dt}\\right|} \\\ '
output+='&\\textrm{Binormal:} &&\\vec{B} = \\vec{T}\\times\\vec{N} \\\ '
output+='&\\textrm{Curvature:} &&\\kappa = \\left|\\frac{d\\vec{T}}{ds}\\right| = \\frac{\\left|\\frac{d\\vec{T}}{dt}\\right|}{v} = \\frac{|\\vec{v}\\times\\vec{a}|}{v^3} \\\ '
output+='&\\textrm{Torsion:} &&\\tau = \\frac{(\\vec{v}\\times\\vec{a})\\cdot \\frac{d^3\\vec{p}}{dt^3}}{|\\vec{v}\\times\\vec{a}|^2} \\\ '
output+='&\\textrm{Frenet-Serret formulas:} &&\\begin{cases} \\frac{d\\vec{T}}{ds} = \\kappa\\vec{N} \\\ \\frac{d\\vec{N}}{ds} = -\\kappa\\vec{T} + \\tau\\vec{B} \\\ \\frac{d\\vec{B}}{ds} = -\\tau\\vec{N} \\end{cases} \\\ '
output+='&\\textrm{Osculating circle:} &&\\textrm{Circle of radius}\\; 1/\\kappa \\; \\textrm{in the plane through}\\;\\vec{p}\\;\\textrm{with normal vector}\\;\\vec{B}\\textrm{, centered at}\\; \\vec{p} + \\bigl(1/\\kappa\\bigr)\\vec{N}. '
output+=' \\end{aligned} $$'

html(output)


@interact
def Motion(curve=input_box(default='[(1/15)*(t^2+1)*cos(t),(1/15)*(t^2+1)*sin(t),t]',type=str,width=100,label='$[x(t),y(t),z(t)]=$'),\
           interval=input_box(default='[0,3*pi]',type=str,width=30,label='$[t_{\\min},t_{\\max}]=$'),\
           param=slider(0,100,1,label='Parameter $t_0=$',display_value=False),\
           show_position=('Show position (blue point):',True),\
           show_velocity=('Show velocity (green):',False),\
           show_acceleration=('Show acceleration (pink) with tangential (green) and normal (red) projections:',False),\
           show_TNBframe=('Show TNB-frame (green/red/blue):',False),\
           show_osccircle=('Show osculating circle:',False),\
           show_computations=('Show computations:',True)):

    [tmin,tmax]=sage_eval(interval)
    ndigits=5 # Numerical precision up to 5 digits behind the decimal dot

    #
    # Vector and scalar functions pertaining to the motion
    #
    position(t)=sage_eval(curve,locals={'t':t})
    velocity = derivative(position,t)
    acceleration = derivative(velocity,t)
    jerk = derivative(acceleration,t)
    speed = sqrt(velocity[0]^2+velocity[1]^2+velocity[2]^2)
    tangent = (1/speed)*velocity
    dT = derivative(tangent,t)
    dTnorm = sqrt(dT[0]^2+dT[1]^2+dT[2]^2)
    normal = (1/(dTnorm))*dT
    binormal = tangent.cross_product(normal)
    vcrossa=velocity.cross_product(acceleration)
    vcrossanorm=sqrt(vcrossa[0]^2+vcrossa[1]^2+vcrossa[2]^2)
    curvaturefunction = (1/speed^3)*vcrossanorm
    torsionfunction = (1/(vcrossanorm^2))*(vcrossa.dot_product(jerk))

    picture=parametric_plot3d(position(t),(t,tmin,tmax),color='red',thickness=3)
    plotsize=(vector(picture.bounding_box()[1])-vector(picture.bounding_box()[0])).norm()

    #
    # Scalar quantities at point
    #
    t0=RR(tmin+(param/100)*(tmax-tmin))
    pos_tzero = position(t0)
    tangent_comp = (1/speed(t0))*velocity(t0).dot_product(acceleration(t0))
    normal_comp = (1/speed(t0))*vcrossanorm(t0)
    curvature = curvaturefunction(t0)
    torsion = torsionfunction(t0)

    #
    # Vector, plane, and circle plot
    #
    pos = point3d(pos_tzero, color='blue',size=10)
    vel = arrow3d(pos_tzero, pos_tzero + velocity(t0), rgbcolor=(0,0.5,0))
    tan = arrow3d(pos_tzero, pos_tzero + tangent(t0), rgbcolor=(0,1,0) )
    nor = arrow3d(pos_tzero, pos_tzero + normal(t0), rgbcolor=(0.5,0,0))
    bin = arrow3d(pos_tzero, pos_tzero + binormal(t0), rgbcolor=(0,0,0.5))
    acc = arrow3d(pos_tzero, pos_tzero + acceleration(t0), rgbcolor=(1,0,1))
    tanproj = arrow3d(pos_tzero, pos_tzero + tangent_comp*tangent(t0), rgbcolor=(0,1,0))
    norproj = arrow3d(pos_tzero, pos_tzero + normal_comp*normal(t0), rgbcolor=(0.5,0,0))
    TNplane = parametric_plot3d(pos_tzero+u*tangent_comp*tangent(t0)+v*normal_comp*normal(t0),\
                                (u,0,1),(v,0,1),opacity=0.2,color='yellow')
    oscdisk = parametric_plot3d(pos_tzero+(1/curvature)*(1+v*sin(u))*normal(t0)\
        +(1/curvature)*v*cos(u)*tangent(t0),(u,0,2*pi),(v,0,1),color='lightblue',thickness=3,opacity=0.5)
    osccircle = parametric_plot3d(pos_tzero+(1/curvature)*(1+sin(u))*normal(t0)\
        +(1/curvature)*cos(u)*tangent(t0),(u,0,2*pi),color='darkblue',thickness=3,opacity=0.8)

    #
    # Assemble plot
    #
    if (show_position == True):
        picture+=pos
    if (show_velocity == True):
        picture+=vel
    if (show_acceleration == True):
        picture+=acc+tanproj+norproj+TNplane
    if (show_TNBframe == True):
        picture+=tan+nor+bin
    if (show_osccircle == True):
        picture+=oscdisk+osccircle

    #
    # Assemble written output
    #
    output='$$ \\begin{aligned} '
    output+='&\\textrm{Position function:} &&\\vec{p}(t) = \\left\\langle %s,%s,%s \\right\\rangle \\quad \\textrm{for}\\; t \\in [%s,%s] \\\ '%(latex(position[0](t)),latex(position[1](t)),latex(position[2](t)),latex(tmin),latex(tmax))
    output+='&\\textrm{Parameter value:} &&t_0 = %s \\\ '%(latex(t0.n(digits=ndigits)))
    output+='&\\textrm{Position:} &&\\vec{p}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle \\\ '%(latex(pos_tzero[0].n(digits=ndigits)),latex(pos_tzero[1].n(digits=ndigits)),latex(pos_tzero[2].n(digits=ndigits)))
    output+='&\\textrm{Velocity:} &&\\vec{v}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle \\\ '%(latex(velocity[0](t0).n(digits=ndigits)),latex(velocity[1](t0).n(digits=ndigits)),latex(velocity[2](t0).n(digits=ndigits)))
    output+='&\\textrm{Speed:} &&v(t_0) = %s \\\ '%(latex(speed(t0).n(digits=ndigits)))
    output+='&\\textrm{Acceleration:} &&\\vec{a}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle = %s \\cdot\\vec{T}(t_0) + %s \\cdot\\vec{N}(t_0) \\\ '%(latex(acceleration[0](t0).n(digits=ndigits)),latex(acceleration[1](t0).n(digits=ndigits)),latex(acceleration[2](t0).n(digits=ndigits)),latex(tangent_comp.n(digits=ndigits)),latex(normal_comp.n(digits=ndigits)))
    output+='&\\textrm{Unit tangent:} &&\\vec{T}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle \\\ '%(latex(tangent[0](t0).n(digits=ndigits)),latex(tangent[1](t0).n(digits=ndigits)),latex(tangent[2](t0).n(digits=ndigits)))
    output+='&\\textrm{Unit normal:} &&\\vec{N}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle \\\ '%(latex(normal[0](t0).n(digits=ndigits)),latex(normal[1](t0).n(digits=ndigits)),latex(normal[2](t0).n(digits=ndigits)))
    output+='&\\textrm{Binormal:} &&\\vec{B}(t_0) = \\left\\langle %s,%s,%s \\right\\rangle \\\ '%(latex(binormal[0](t0).n(digits=ndigits)),latex(binormal[1](t0).n(digits=ndigits)),latex(binormal[2](t0).n(digits=ndigits)))
    output+='&\\textrm{Curvature:} &&\\kappa(t_0) = %s \\\ '%(latex(curvature.n(digits=ndigits)))
    output+='&\\textrm{Torsion:} &&\\tau(t_0) = %s '%(latex(torsion.n(digits=ndigits)))
    output+=' \\end{aligned} $$'

    #
    # Show written output and picture
    #
    if (show_computations == True):
        html(output)
    show(picture,aspect_ratio=[1,1,1])

Motion in 3-dimensional space

$$\begin{aligned} &\textrm{Position:} &&\vec{p}(t) = \langle x(t),y(t),z(t) \rangle \quad \textrm{for } t \in [t_{\min},t_{\max}] \\ &\textrm{Velocity:} &&\vec{v} = \frac{d\vec{p}}{dt} \\ &\textrm{Speed:} &&v = |\vec{v}| \\ &\textrm{Acceleration:} &&\vec{a} = \frac{d\vec{v}}{dt} = \frac{d^2\vec{p}}{dt^2} = a_{T}\vec{T} + a_N\vec{N} = \left(\frac{dv}{dt}\right) \vec{T} + \left(\kappa v^2\right) \vec{N} = \left(\frac{\vec{v}\cdot\vec{a}}{v}\right)\vec{T} + \left(\frac{|\vec{v}\times\vec{a}|}{v}\right)\vec{N} \\ &\textrm{Unit tangent:} &&\vec{T} = \frac{\vec{v}}{v} \\ &\textrm{Unit normal:} &&\vec{N} = \frac{\frac{d\vec{T}}{dt}}{\left|\frac{d\vec{T}}{dt}\right|} \\ &\textrm{Binormal:} &&\vec{B} = \vec{T}\times\vec{N} \\ &\textrm{Curvature:} &&\kappa = \left|\frac{d\vec{T}}{ds}\right| = \frac{\left|\frac{d\vec{T}}{dt}\right|}{v} = \frac{|\vec{v}\times\vec{a}|}{v^3} \\ &\textrm{Torsion:} &&\tau = \frac{(\vec{v}\times\vec{a})\cdot \frac{d^3\vec{p}}{dt^3}}{|\vec{v}\times\vec{a}|^2} \\ &\textrm{Frenet-Serret formulas:} &&\begin{cases} \frac{d\vec{T}}{ds} = \kappa\vec{N} \\ \frac{d\vec{N}}{ds} = -\kappa\vec{T} + \tau\vec{B} \\ \frac{d\vec{B}}{ds} = -\tau\vec{N} \end{cases} \\ &\textrm{Osculating circle:} &&\textrm{Circle of radius}\; 1/\kappa \; \textrm{in the plane through}\;\vec{p}\;\textrm{with normal vector}\;\vec{B}\textrm{, centered at}\; \vec{p} + \bigl(1/\kappa\bigr)\vec{N}. \end{aligned}$$

Interact: please open in CoCalc