Making an iterative spring compared to this non-iterative spring?

So I made this non-iterative spring (spring that uses time to calculate what spring movements will look like, i think thats what that means anyway, and I wanted to know how you would do this iteratively instead. By this, I mean you won’t need the time argument only the spring variables themself.

My current non-iterative code (well, using time):

--Physical parts, to see the movement of the spring
local target = workspace.Target
local current = workspace.Current

local function create(c, A, B)
	--New table for the spring
	local s = {}
	
	--A for amplitude, B for frequency
	s.A = A
	s.B = B
	
	--c is the starting value, t is the target value, and s is the current value
	s.c = c
	s.t = c
	s.s = c
	
	--Function which updates the spring, where x is the time elapsed since the target was updated
	s.update = function(x)
		s.s = (s.c - s.t) * math.cos(s.B * x) * math.exp(-s.A * x) + s.t
	end
	return s
end

--Create the spring, and set time to 0
local spring = create(0, 4, 5)
local d = 0

--Update spring every frame
game:GetService("RunService").RenderStepped:Connect(function(s)
	--d is time elapsed
	d = d + s
	spring.update(d)
	--set physical positions
	target.Position = Vector3.new(spring.t, 0, 1.05)
	current.Position = Vector3.new(spring.s, 0, 1.1)
end)

--Wait before updating the spring first
wait(2)

d = 0
spring.c = spring.s
spring.t = 10

wait(2)

d = 0
spring.c = spring.s
spring.t = 0

Wondering about how to go about this iteratively instead, would like some help on this.

The non-iterative method is much easier, because you don’t have to find the integrals of the spring force to obtain the velocity and position. I’ll attempt to do so here:

Assuming that Hooke’s law holds (it does except under extreme extension / compression) then the force exerted by a spring is:

F = k * Δx

Where k is the spring constant (different materials and shape affect it) and Δx is the offset from the resting position. Integrating this is simple, but instead of integrating relative to the offset, we need to integrate relative to time, t. This gives us this equation:

V = k * Δx * t + C

C in this case would be the previous calculated velocity, so we can write it like so:

V = v0 + k * Δx * t

When we integrate this to get the position, we have:

P = v0 * t + (k * Δx * t^2)/2 + C

And C in this case would be the initial spring height, p0:

P = p0 + v0 * t + (k * Δx * t^2)/2

So, using p0 and v0 from the previous iteration we can find the current Δx (difference between p0 and the resting height) and use the time since the last iteration (t) to find the new spring velocity and position.

local restingHeight = 0
local springConstant = 0.5
local function spring(lastHeight, lastVelocity, deltaTime)
    -- When above, we want the force to be negative to bring it back down
    local heightOffset = restingHeight - lastHeight
    local deltaVelocity = springConstant * heightOffset * deltaTime
    local height = lastHeight + (lastVelocity + deltaVelocity / 2) * deltaTime
    return height, lastVelocity + deltaVelocity
end

Edit: I’m back and had some time to write an example.

I don’t think this is correct, because kΔx is dependent on t, so your integration is wrong. Also, OP appears to be using some sort of damped harmonic oscillator model (note the exponential decay), so Hooke’s law doesn’t work here.

@OP: Since you already have the analytical solution, the easiest method to make an iterative version would be to calculate the derivative of the non-iterative function. The derivative multiplied by the input increment (time increment in this case) equals the output increment (change in position).

Also, note that amplitude isn’t the right word for your ‘A’ parameter. You should call it the decay rate. The amplitude would be s.t - s.c.

I’m assuming you want to switch from using an explicit formula to a recursive one so that you can eliminate time. If so, you have to use (up to) the second derivative as that yields a linear homogeneous differential equation (no time variable in the equation). This will be the original force formula that was used to derive your position time function.

I think your original force equation used to derive that must of had linear dampening where the result was an underdamped spring. The force equation (assume mass = 1) is:

F = -k*(x - x0) - m*velocity

k is the spring coefficient, and, m is the drag coefficient.

From there, using semi-implicit euler as the integrator (pretty good but can drift out of phase),

velocity += acceleration * dt
position += velocity * dt

so

velocity = velocity + (-k*(position - restingPosition) - m*velocity)*dt
position = position + velocity*dt

dt = 1/100 is pretty accurate, however, if you want more accuracy, I would suggest using verlet as the integration scheme.

EDIT:

From your formula,
the initial position (before you begin iteration) should be position = c
the initial velocity (before you begin iteration) should be velocity = 0 (starts from rest)
restingPosition = t (target value)
position = s (current value)

Also, if you want to set amplitude and frequency directly (as done in your explicit formula) then,
Using A for amplitude, B for frequency

m = 2*A

and

k = A^2 + B^2

The force exerted isn’t dependent on time, only the extension/compression of the spring. The units for Force is newtons (N) which is m/s^2. Since t is squared in the position function, it gives us the position in meters (studs in our case).

Ah, yeah, that does put a dampener on this method xD

I’ll transform this into

a * cos(f * x) / e^x + h

I moved the decay rate out of the exponential decay function to become a part of the amplitude. The new a = Δh/e^(s.A) if you want to convert a rate of decay in the old model to the amplitude in this model. Making that substitution doesn’t change the function at all. The function that I gave equates to this portion of the equation:

a * cos(f * x)

Where a and f are hidden in the definition of the difference from resting height Δx and spring constant k. The wave effect of cosine is hidden in Δx as well because as the velocity carries the height below equilibrium, this value becomes negative. I tested this by running this code:

local restingHeight = 0
local springConstant = 1
local function spring(lastHeight, lastVelocity, deltaTime)
    -- When above, we want the force to be negative to bring it back down
    local heightOffset = restingHeight - lastHeight
    local deltaVelocity = (springConstant * heightOffset) * deltaTime
    local height = lastHeight + (lastVelocity + deltaVelocity / 2) * deltaTime
    return height, lastVelocity + deltaVelocity
end

local lastHeight = 1
local lastVelocity = 0
local stepsPerFrame = 4
local framesPerSec = 60
local stepsPerSec = framesPerSec * stepsPerFrame
local seconds = 10
for i = 1, stepsPerSec * seconds do
    print("Height:   " .. lastHeight)
    lastHeight, lastVelocity = spring(lastHeight, lastVelocity, 1/stepsPerSec)
end
print("Height:   " .. lastHeight)
print("Velocity: " ..lastVelocity)

What I found was that as the number of steps per a second increased, the valleys and crests seem to near -1 and 1, indeed it seems to converge to a sinusoidal wave. I also noted that if the step was not small enough, then the spring would actually diverge. I picked 4 steps a second and 60 frames a second to emulate what the Roblox physics engine uses. 1 step every 60 frames is what you would get by connecting to the heartbeat event, and that works just fine.

The function I gave you however is not damped by 1/e^x and will continue to bounce up and down. To add a linear dampening to the equation (such as friction) we simply make the force:

F = k * Δx - d * v0

Where v0 is the last velocity. This is because the force of friction is linearly dependent on velocity. Following this all the way down, we get these equations for the velocity and position after simplifying:

V = v0 + (k * Δx - d * v0) * t

and

P = p0 + v0 * t + (k * Δx - d * v0) / 2 * t^2

I threw this into my code above and got a nice damped spring. You will have to play with the spring constant to manipulate the period of the spring and the dampener constant to manipulate how long it takes for the spring to reach near equilibrium. I’m not sure what value make the dampening critical.

Here is the end result that you can see for yourself:

local restingHeight = 0
local springConstant = 1
local dampening = 0.1
local function spring(lastHeight, lastVelocity, deltaTime)
    -- When above, we want the force to be negative to bring it back down
    local heightOffset = restingHeight - lastHeight
    local deltaVelocity = (springConstant * heightOffset - dampening * lastVelocity) * deltaTime
    local height = lastHeight + (lastVelocity + deltaVelocity / 2) * deltaTime
    return height, lastVelocity + deltaVelocity
end

local lastHeight = 1
local lastVelocity = 0
local stepsPerFrame = 4
local framesPerSec = 60
local stepsPerSec = framesPerSec * stepsPerFrame
local seconds = 10
for i = 1, stepsPerSec * seconds do
    print("Height:   " .. lastHeight)
    lastHeight, lastVelocity = spring(lastHeight, lastVelocity, 1/stepsPerSec)
end
print("Height:   " .. lastHeight)
print("Velocity: " ..lastVelocity)

Edit: v0 and v were the same value in the equations above. I’ve made the change.

1 Like

The extension and compression are dependent on time, because with force being mass times acceleration, the force induces a change in the extension/compression of the spring, which means that the extension/compression is changing over time, and therefore the force also changes. The force is not a static property of the spring — it is a component of the spring’s state which evolves over time in parallel with other components of the state, namely displacement and velocity.

Ah, I see what you are saying now. My equation for the force factors this into Δx, the difference between the current position and the resting position (in essence the extension / compression). Hookes law says that a spring’s force is linearly dependent on this value by a constant except for when extremely compressed or extended. You can see my equation for the force of a spring on the Wikipedia page on Hooke’s law under “For Linear Springs”.

If you understand what I’m saying, then allow me to state that as a direct consequence of Δx‘s dependence on time, the integral of kΔx is not simply kΔx * t + C, although this will give a good approximation in a small neighborhood around t=0. You cannot solve Hooke’s law by direct integration. On the contrary, if Δx = cos(t), then F = -m cos(t), and mv = ʃF dt = m dx/dt = -m sin(t), which is different than F * t. In general, d(Ft)/dt = F + t dF/dt.

Since your equations are accurate for values in a small neighborhood of t=0, it turns out that they are usable for iterative solutions. This is why your code goes haywire when the time increment becomes larger, because the higher-order terms in the Taylor series become relevant. This is the nature of iterative solutions, so this isn’t necessarily a bad thing. However, it should be said that your equations are incorrect as stated.

2 Likes