Minimum distance travelled in order to reach a known maximum acceleration given a constant jerk value

calculusphysics

Problem

Given the following:

$A$ = maximum acceleration.
$J$ = constant jerk (change in acceleration with respect to time).
Initial velocity and acceleration of 0.

How might one determine the minimum distance that must be travelled in order to reach the maximum acceleration?

Attempt

In the same way that one might solve distance $D$ given a maximum velocity $V$ and constant acceleration $A$ like so:

$D = \frac{V^2}{2A}$

I thought I might be able to apply the same to my problem by swapping out maximum velocity with maximum acceleration, and the constant acceleration for constant jerk, before realising that this seems to give the velocity at which the maximum acceleration is reached, and not the distance.

$V = \frac{A^2}{2J}$

I'm as yet unsure how to apply what I know so far to determine the distance travelled before the maximum acceleration is reached. Any help would be greatly appreciated!

Application Context

I'm currently attempting to improve a linear motion profile for some linear actuators that are present in an art installation that I'm working on in a small team.

Previously, the linear actuator was driven using a constant acceleration selected by us. We could achieve fast movement nicely with high acceleration, and we could achieve slow movement nicely with slow acceleration, however slow movement with a high acceleration resulted in a lot of "chatter" from the motor.

The idea is to use a constant jerk rather than a constant acceleration to achieve smoother behaviour at both low and high speeds. The implementation seems to be almost complete, however I'm struggling to determine the distance at which I should begin decreasing deceleration in order to arrive at the target position with a velocity and acceleration of 0. It occurs to me that this is the same as solving the distance required to reach the maximum acceleration, hence the framing of this question.

Best Answer

If $\dot{a}=J$, then repeated integration yields

$$ a(t)=Jt+C_1 \xrightarrow{a(0)=0} a(t)=Jt$$ $$ v(t)=\frac{Jt^2}{2}+C_2 \xrightarrow{v(0)=0} v(t)=\frac{Jt^2}{2}$$ $$ x(t)=\frac{Jt^3}{6}+C_3 \xrightarrow{x(0)=0}x(t)=\frac{Jt^3}{6}$$

Solving the last equation for $t$

$$ t=\left(\frac{6x}{J}\right)^{1/3}$$

Substitute in $a(t)$

$$ a=J^{2/3}(6x)^{1/3}$$

Solve for $x$

$$ x(a)=\frac{a^3}{6J^2}$$

Related Solutions

Calculate the constant acceleration needed for the discrete time trajectory to intersect a given target point

A quick background refresher:

In continuous time, the velocity $v(t)$ and position $x(t)$ with a constant acceleration $a$ are given by $$\left\lbrace\begin{aligned} v(t) &= \int_{0}^{t} a \, d \tau = v_0 + a t \\ x(t) &= \int_{0}^{t} v(\tau) \, d \tau = x_0 + v_0 t + \frac{1}{2} a t^2 \end{aligned}\right.\tag{1}\label{NA1}$$

In discrete time, the velocity and position are given by difference equations. There are three most common definitions:

A. Position updated before velocity: $$\left\lbrace\begin{aligned} x_{i+1} &= x_{i} + v_{i} \\ v_{i+1} &= v_{i} + a \\ \end{aligned}\right. \implies x_N = x_0 + N v_0 + \frac{N (N - 1)}{2} a \tag{2a}\label{NA2a}$$

B. Velocity updated before position: $$\left\lbrace\begin{aligned} v_{i+1} &= v_{i} + a \\ x_{i+1} &= x_{i} + v_{i+1} \\ \end{aligned}\right. \implies x_N = x_0 + N v_0 + \frac{N(N+1)}{2} a \tag{2b}\label{NA2b}$$

C. Velocity updated before and after position: $$\left\lbrace\begin{aligned} v_{i+\frac{1}{2}} &= v_{i} + \frac{1}{2} a \\ x_{i+1} &= x_{i} + v_{i+\frac{1}{2}} \\ v_{i+1} &= v_{i+\frac{1}{2}} + \frac{1}{2} a \\ \end{aligned}\right. \implies x_N = x_0 + N v_0 + \frac{N^2}{2} a \tag{2c}\label{NA2c}$$

In this last case $\eqref{NA2c}$, velocity updates are essentially half a time step out of sync with the position updates, but the position function is most similar to the continuous time one, $\eqref{NA1}$.

(While the formulae above are shown in one dimension, you can use them in 2D or 3D just as well. The time steps ($i$ and $N$) are the same (shared) across all dimensions, but each dimension has their own acceleration $a$. The actual magnitude of the acceleration is then $\sqrt{a_x^2 + a_y^2}$ for 2D, and $\sqrt{a_x^2 + a_y^2 + a_z^2}$ for 3D.)

The problem at hand is to find the smallest $0 \le N \in \mathbb{N}$ and $a_x$ and $a_y$, when $x_0$, $y_0$, $v_{x0}$, $v_{y0}$ are known, but acceleration is limited, $a_x^2 + a_y^2 \le A_{max}^2$.

If we solve $\eqref{NA2a}$ for $a$ for both dimensions, we get $$\left\lbrace\begin{aligned} a_x &= \frac{2 (x_N - x_0 - N v_{x0})}{N (N - 1)} \\ a_y &= \frac{2 (y_N - y_0 - N v_{y0})}{N (N - 1)} \\ \end{aligned}\right . \tag{3a}\label{NA3a}$$ Solving $\eqref{NA2b}$ similarly yields $$\left\lbrace\begin{aligned} a_x &= \frac{2 (x_N - x_0 - N v_{x0})}{N (N + 1)} \\ a_y &= \frac{2 (y_N - y_0 - N v_{y0})}{N (N + 1)} \\ \end{aligned}\right . \tag{3b}\label{NA3b}$$ and $\eqref{NA2c}$ yields $$\left\lbrace\begin{aligned} a_x &= \frac{2 (x_N - x_0 - N v_{x0})}{N^2} \\ a_y &= \frac{2 (y_N - y_0 - N v_{y0})}{N^2} \\ \end{aligned}\right . \tag{3c}\label{NA3c}$$ To find out the approximate number of steps (which we need to round upwards to not exceed the specified acceleration $A_{max}$), we could try and solve $$a_x^2 + a_y^2 \le A_{max}^2 \tag{4}\label{NA4}$$ for $N \in \mathbb{R}$ (so the actual number of time steps needed would be rounded up, $\left\lceil N \right\rceil$).

Unfortunately, while there is an algebraic solution (it is a form of a fourth-degree function that does have algebraic solutions), it is too complicated to be useful.

Fortunately, we can use the bisection method (or binary search) to very efficiently find the smallest $N$, and therefore also the largest acceleration $(a_x , a_y)$ not exceeding $A_{max}$ in magnitude, that reaches the target coordinates in minimum number of time steps.

Note that in practice, square root and division operations take much more time than addition or multiplication. Because the bisection method/binary search will have to evaluate the same expression potentially dozens of times, it is a good idea to "optimize" $\eqref{NA4}$ first.

If we use $\Delta_x = x_N - x_0$ and $\Delta_y = y_N - y_0$, then $\eqref{NA4}$ optimizes to: $$(\Delta_x - N v_{x0})^2 + (\Delta_y - N v_{y0})^2 - N^2 (N-1)^2 \frac{A_{max}^2}{4} \le 0 \tag{5a}\label{NA5a}$$ $$(\Delta_x - N v_{x0})^2 + (\Delta_y - N v_{y0})^2 - N^2 (N+1)^2 \frac{A_{max}^2}{4} \le 0 \tag{5b}\label{NA5b}$$ $$(\Delta_x - N v_{x0})^2 + (\Delta_y - N v_{y0})^2 - N^4 \frac{A_{max}^2}{4} \le 0 \tag{5c}\label{NA5c}$$ depending on how the velocity is updated with respect to position updates, respectively. (This means each iteration when finding $N$ takes only something like 8 multiplications and 5 additions or subtractions. Since the number of iterations needed is roughly $2 \log_2 N$, this is actually pretty efficient. Even if we find $N \approx 10,000$, we only end up doing a couple of hundred multiplications and a hundred and thirty additions or subtractions. The algebraic solution has more terms than that!)

Furthermore, the optimized form is defined even for $N = 0$ and $N = 1$, so there is no risk of accidental divide-by-zero error. This makes the bisection search easier to write, too.

If we use EQ5(N) for $\eqref{NA5a}$/$\eqref{NA5b}$/$\eqref{NA5c}$, the pseudocode to find $N$ is something like

Function StepsNeeded(x0,y0, xN,yN, vx,vy, Amax):
    Let  dx = xN - x0
    Let  dy = yN - y0
    Let  Nmin = 0
    Let  Nmax = 2

    # Find the Nmin,Nmax range spanning the solution
    While EQ5(Nmax) > 0:
        Let  Nmin = Nmax
        Let  Nmax = Nmax * 2
    End While

    Let  N = Nmax
    While Nmax > Nmin:
        Let  N = Nmax - floor((Nmax - Nmin)/2)
        Let  C = EQ5(N)
        If C > 0:
            Let  Nmin = N
        Else If C < 0:
            If Nmax > N:
                Let  Nmax = N
            Else:
                Break While
            End If
        Else:
            Break loop
        End If
    End While

    Return N
End Function

When we do have $N$, it is just a matter of plugging it back to $\eqref{NA3a}$/$\eqref{NA3b}$/$\eqref{NA3c}$ to obtain the acceleration vector $(a_x , a_y)$.

The first While loop finds the smallest range containing a root, so that Nmax is a sufficient number of steps (yields an acceleration that does not exceed the given limit), and Nmin is an insufficient number of steps (yields an acceleration that exceeds the given limit). For each iteration $1 \le i \in \mathbb{N}$, Nmin= $2^{i-1}$ (except Nmin=0 for $i = 1$, and Nmin=2 for $i = 2$) and Nmax=$2^i$, so only a few iterations are ever needed.

The second While loop uses the bisection method (or binary search) to find the root within that range. As it halves the range on each iteration, it does at most as many iterations as the first loop did.

Because we want an $N$ that yields $a_x^2 + a_y^2 \le A_{max}^2$, we exclude the minimum (Nmin) and include the maximum (Nmax) in the search. This is why we use N = Nmax - floor((Nmax - Nmin)/2), too: it never yields Nmin. This is also why we start with Nmin = 0: the smallest N returned is 1.

Here is an example Python (works in both Python 2 and Python 3) program, that solves the number of steps and the acceleration, when given the starting x and y coordinates, target x and y coordinates, initial velocity x and y components, and the maximum acceleration allowed:

from math import sqrt
from sys import stderr, argv, exit

# This file is in public domain. No guarantees, no warranties.
# Written by Nominal Animal <[email protected]>.

def solveacceleration(start, finish, velocity, maxaccel):
    dx = finish[0] - start[0]
    dy = finish[1] - start[1]
    vx = velocity[0]
    vy = velocity[1]
    ac = 0.25 * maxaccel * maxaccel

    if maxaccel <= 0.0:
        return 0, (0.0, 0.0)

    nmin = 0
    nmax = 2
    while True:
        cx = dx - nmax*vx
        cy = dy - nmax*vy

        # A) cc = cx**2 + cy**2 - ac * nmax**2 * (nmax-1)**2
        # B) cc = cx**2 + cy**2 - ac * nmax**2 * (nmax+1)**2
        # C) cc = cx**2 + cy**2 - ac * nmax**4

        cc = cx**2 + cy**2 - ac * nmax**2 * (nmax-1)**2
        # stderr.write("nmin=%u, nmax=%u, cc=%.6f\n" % (nmin, nmax, cc))
        if cc > 0:
            nmin = nmax
            nmax = nmax * 2
        else:
            break

    n = nmax
    while nmax > nmin:
        n = nmax - int((nmax - nmin) / 2)
        cx = dx - n*vx
        cy = dy - n*vy

        # A) cc = cx**2 + cy**2 - ac * n**2 * (n-1)**2
        # B) cc = cx**2 + cy**2 - ac * n**2 * (n+1)**2
        # C) cc = cx**2 + cy**2 - ac * n**4
        cc = cx**2 + cy**2 - ac * n**2 * (n-1)**2

        # stderr.write("nmin=%u, nmax=%u, n=%u, cc=%.6f\n" % (nmin, nmax, n, cc))
        if cc > 0:
            nmin = n
        elif cc < 0:
            if nmax == n:
                break
            else:
                nmax = n
        else:
            break

    # A) a_ = 2*(d_ - n*v_) / (n * (n - 1))
    # B) a_ = 2*(d_ - n*v_) / (n * (n + 1))
    # C) a_ = 2*(d_ - n*v_) / (n**2)
    ax = 2*(dx - n*vx) / (n * (n - 1))
    ay = 2*(dy - n*vy) / (n * (n - 1))
    return n, (ax, ay)


if __name__ == '__main__':

    if len(argv) != 8:
        stderr.write("\n")
        stderr.write("Usage: python %s x0 y0 x1 y1 xv yv amax\n" % argv[0])
        stderr.write("\n")
        exit(0)

    p0 = ( float(argv[1]), float(argv[2]) )
    p1 = ( float(argv[3]), float(argv[4]) )
    v0 = ( float(argv[5]), float(argv[6]) )
    amax = float(argv[7])

    n, a = solveacceleration(p0, p1, v0, amax)

    print("#  From: %9.3f %9.3f" % p0)
    print("#    To: %9.3f %9.3f" % p1)
    print("#    v0: %+9.3f %+9.3f" % v0)
    print("# Steps: %9d" % n)
    print("#     a: %+9.3f %+9.3f  (%.3f of maximum %.3f)" %
          (a[0], a[1], sqrt(a[0]*a[0]+a[1]*a[1]), amax))

    p = p0
    v = v0
    print("# step    x         y          vx        vy")
    print("")
    for i in range(0, n+1):
        # A) print("%-5d %9.3f %9.3f   %+9.3f %+9.3f" % (i, p[0], p[1], v[0], v[1]))
        #    p = (p[0]+v[0], p[1]+v[1])
        #    v = (v[0]+a[0], v[1]+a[1])

        # B) v = (v[0]+a[0], v[1]+a[1])
        #    print("%-5d %9.3f %9.3f   %+9.3f %+9.3f" % (i, p[0], p[1], v[0], v[1]))
        #    p = (p[0]+v[0], p[1]+v[1])

        # C) v = (v[0]+0.5*a[0], v[1]+0.5*a[1])
        #    print("%-5d %9.3f %9.3f   %+9.3f %+9.3f" % (i, p[0], p[1], v[0], v[1]))
        #    p = (p[0]+v[0], p[1]+v[1])
        #    v = (v[0]+0.5*a[0], v[1]+0.5*a[1])

        print("%-5d %9.3f %9.3f   %+9.3f %+9.3f" % (i, p[0], p[1], v[0], v[1]))
        p = (p[0]+v[0], p[1]+v[1])
        v = (v[0]+a[0], v[1]+a[1])

The solveacceleration() function takes the starting point, target or finish point, and the velocity as two-component tuples. If you include the commented out stderr.write() lines, the function will output a line per iteration when looking for the answer; you'll see that it does very few iterations even for the most complex cases. It returns the number of steps, plus the acceleration as a two-component tuple.

The program itself outputs the parameters and the solution, including the location and velocity at each step, in a format suitable for e.g. gnuplot. If you save the output as e.g out.txt, you can use

plot "out.txt" u 2:3 notitle w lines lc -1, \
     "out.txt" u 2:3:1 notitle w labels \
               hypertext point pt 6

in gnuplot to draw the trajectory, with small circles around each time step; hovering over each circle shows you the time step as a tooltip. To include the velocity in the tooltip, use

plot "out.txt" u 2:3 notitle w lines lc -1, \
     "out.txt" u 2:3:(sprintf("%s: (%s,%s)", stringcolumn(1), stringcolumn(4), stringcolumn(5))) \
               notitle w labels hypertext point pt 6

As it is written above, the example uses the $\eqref{NA2a}$ ("A") logic, position updated before velocity. I've commented the sections where you need to modify the code for the other two ("B" or "C").

Time to travel a given distance, given two separate, but constant acceleration and deceleration rates and a max velocity.

You need to figure out where the breakpoints are, then it just becomes solving constant acceleration between the breakpoints. For your first example, you probably accelerate up to maximum speed, stay at that speed until you need to decelerate, and decelerate fully. The only reason this would not be true is if you can't get up to maximum speed before you need to decelerate, so take that assumption. Now see if you can get to maximum speed before the atmosphere starts. If you can, you just accelerate in vacuum up to maximum speed. See if you can decelerate from maximum speed while within the atmosphere. If you can, do that, and you have the breakpoints where you stop accelerating and start decelerating. If you can't decelerate within the atmosphere, decelerate fully in the atmosphere and solve for the velocity you need at the top of the atmosphere. Start decelerating in vacuum to meet this speed.

If you can't get up to max speed and back down again, you need to locate the altitude of maximum velocity and the velocity at that altitude. You will accelerate fully until you hit that point, then decelerate fully. Guess whether it is in or out of the atmosphere. One of the regimes will just have one acceleration rate, the other will have two.

For the second, accelerate fully in the atmosphere until either you hit top speed or you hit the top of the atmosphere. In the first case, just keep going to $45$ km. In the second, keep accelerating in vacuum until maximum speed, then just hold it.