【此文章原作者为huw bowles，原文链接为http://www.huwbowles.com/spring-dynamics-production/，下面使用google翻译并且人工润色了一些。嗯，物理模拟什么的最好玩了。】
介绍
弹簧是一种非常有用且基本的工具，可以表达各种动力。 在游戏和物理模拟中，我广泛地使用着它们。它们非常适合进行lerps 或 transitions (with under/overshoot).。 它们还是移动刚性体的理想工具-阻尼弹簧在数学上等效于PD控制。
在过去的几年中，我在Studio Gobo编写的弹簧动力学代码有了很大进步。 这篇文章记录了这种进步，希望对其他人有用。
弹簧实现
基础——简单的弹簧实现
这是一个由物理法则直接实现的弹簧系统
class Spring
{
// 弹簧目标
Vector3 m_targetPos;
// 动力状态
Vector3 m_pos, m_vel;
// 参数 - spring constant
float m_kp;
// 参数 - damping constant
float m_kd;
void Update( const float dt )
{
// compute spring and damping forces
Vector3 Fspring = m_kp * (m_targetPos - m_pos);
Vector3 Fdamp = m_kd * -m_vel;
// integrate dynamics
m_pos += m_vel * dt;
m_vel += (Fspring + Fdamp) * dt;
}
};
弹簧具有目标位置，然后将弹力计算为弹簧常数乘以相对于当前位置和目标位置的偏移量（代码第16行）。我们不希望弹簧永远（或根本不）摆动，因此我们通过施加与速度成比例的力来实现阻尼系统（代码第17行）。 最后，在第20和21行更新动态状态（pos和vel）
在继续学习有趣的东西之前，有一些注意事项需要先解决。 首先，代码虽然是向量形式的，但和在浮点数上运行一样容易（类可以在数据类型上进行模板化）。 目标位置是弹簧的输入，可以根据需要在每帧中进行更新。 其次，建模时无需特别清晰明了。 改变速度时，虽然可以将力除以质量再加到速度上，但我觉得这样做并不是很有用，而是依靠其他spring参数来获得所需的效果。
设置目标速度
基本的弹簧代码忽视了一个重要的条件–它假定目标速度为0。如果将弹簧连接到移动的目标，则该位置将始终滞后于目标，因为弹簧将会以0速度到达目标。 另一方面，如果您指定了目标速度和位置，则可以让该点将紧密跟踪目标。 我们可以添加更多代码来实现此功能：
class Spring
{
// target for spring, set as desired by user
Vector3 m_targetPos, m_targetVel;
// dynamic state
Vector3 m_pos, m_vel;
// spring param - spring constant
float m_kp;
// spring param - damping constant
float m_kd;
void Update( const float dt )
{
// compute spring and damping forces
Vector3 Fspring = m_kp * (m_targetPos - m_pos);
Vector3 Fdamp = m_kd * (m_targetVel - m_vel);
// integrate dynamics
m_pos += m_vel * dt;
m_vel += (Fspring + Fdamp) * dt;
}
};
请注意，我们在两个力计算中是非常对称的。 如果物体的位置和速度是已知的，就可以非常平滑地过渡。 如果参考系相对于世界移动，则也需要此参数。 如果弹簧在太空飞船内部，则目标速度应包括太空飞船的速度。
改进参数–用振动频率代替弹簧常数
而弹簧常数m_kp可以根据物理上有意义的参数来表示：无阻尼角频率（Undamped Angular Frequency，UAF）。 该值具有很重要的含义-它是弹簧振动的频率（每秒的跳动次数）。 这是任何人都可以理解的值，并且比任何常量更好。 它的用法如下：
class Spring
{
// target for spring, set as desired by user
Vector3 m_targetPos, m_targetVel;
// dynamic state
Vector3 m_pos, m_vel;
// spring param - Undamped Angular Frequency (UAF) - oscillation frequency
float m_UAF;
// spring param - damping constant
float m_kd;
void Update( const float dt )
{
// compute spring and damping forces
Vector3 Fspring = m_UAF * m_UAF * (m_targetPos - m_pos);
Vector3 Fdamp = m_kd * (m_targetVel - m_vel);
// integrate dynamics
m_pos += m_vel * dt;
m_vel += (Fspring + Fdamp) * dt;
}
};
改进参数–标准化阻尼参数
现在我们将注意力转向第二个参数-阻尼常数m_kd。 每次弹簧强度/频率UAF发生变化时，都需要重新平衡阻尼项以维持所需的下冲/过冲（undershoot/overshoot），这对于可用性/工作流程是个坏消息。 幸运的是，可以通过将参数更改为阻尼比（Damping Ratio，DR）轻松解决此问题。 它具有很好的特性——如果值<1将会让underdamped的弹簧overshoot，而值> 1则会让overdamped的将缓慢到达目标，而临界阻尼值恰好为1，目标将在此处达到极限。 尽可能快地达到目标而不会超过。
class Spring
{
// target for spring, set as desired by user
Vector3 m_targetPos, m_targetVel;
// dynamic state
Vector3 m_pos, m_vel;
// spring param - Undamped Angular Frequency (UAF) - oscillation frequency
float m_UAF;
// spring param - Damping Ratio (DR)
float m_DR;
void Update( const float dt )
{
// compute spring and damping forces
Vector3 Fspring = m_UAF * m_UAF * (m_targetPos - m_pos);
Vector3 Fdamp = 2.0f * m_DR * m_UAF * (m_targetVel - m_vel);
// integrate dynamics
m_pos += m_vel * dt;
m_vel += (Fspring + Fdamp) * dt;
}
};
现在，设计人员/艺术家/程序员只需确定他们希望弹簧运动有多少overshoot/undershoot，并且可以设置此值。
提高稳定性–集成方案
上面我们已经看到了如何计算弹簧力和阻尼力。 这些力越大，误差将越大。 这往往会使弹簧轨迹的峰值过冲，从而为系统增加能量。 如果弹簧力或阻尼力足够大，则可以完全阻止弹簧收敛并导致模拟爆炸。 减少此错误的一种方法是减少时间步长，从而以更多的计算为代价提供更高的准确性，我们将在下一部分中对此进行探讨。 减少错误的另一种方法是更改集成方案，我们在本节中比较一些不同的方案。The diagrams are taken from an Excel spreadsheetattachedto this post.
Euler method
到目前为止，更新代码已实现了Euler方法。 它仅使用当前状态将状态及时向前移动，称为显式更新。 它很容易实现，但与分析轨迹（浅灰色）相比，往往会严重地超出可见范围：
在上述情况下，阻尼比很低，能量隐藏在系统中，并且于仿真有所不同。 稍微提高阻尼比可以解决这个问题，但是如果我们将阻尼比设置得足够高，则会遇到更严重的问题-尖峰振荡，导致系统爆炸：
半隐式欧拉
一个简单的改进是切换pos和vel的更新顺序：
m_vel += (Fspring + Fdamp) * dt;
m_pos += m_vel * dt;
正如维基百科上记录的那样，该方案非常好，我们在低阻尼值下的问题消失了（橙色）：
不幸的是，该方案在阻尼时仍然容易产生剧烈的振荡：
隐式更新
如上所述，阻尼弹簧与PD控制密切相关，这与论文上的结果一样。 Tan et al引入了稳定的PD控制，对于stiff动力学或较大的时间步长是稳定的[1]。 这个想法是根据下一时间步的状态而不是当前状态来计算力，当前状态是使用一阶泰勒级数预测的。 解决了加速问题后，可以归结为中等数量的附加代码：
class Spring
{
// target for spring, set as desired by user
Vector3 m_targetPos, m_targetVel;
// dynamic state
Vector3 m_pos, m_vel;
// spring param - Undamped Angular Frequency (UAF) - oscillation frequency
float m_UAF;
// spring param - Damping Ratio (DR)
float m_DR;
void Update( const float dt )
{
// compute spring and damping forces
Vector3 posNextState = m_pos + dt * m_vel;
Vector3 Fspring = m_UAF * m_UAF * (m_targetPos - posNextState);
Vector3 Fdamp = 2.0f * m_DR * m_UAF * (m_targetVel - m_vel);
// integrate dynamics
m_pos += m_vel * dt;
m_vel += (Fspring + Fdamp) * dt / (1.0f + dt * 2.0f * m_DR * m_UAF);
}
};
对于低阻尼值，我们可以节约一点，而解决方案（绿色）似乎更符合地面实况（灰色）的频率：
而且，在过阻尼的情况下，我们也会表现出良好的行为：
这些集成方案已实现到随附的电子表格中。 它们也在此处的交互式ShaderToy中运行[2]。
提高顽健性–固定或修复时间步
上面的隐式更新阻止了动态爆炸。 但是，它不能保证弹簧的性能得到保留。 可以看到，在模拟的（绿色）和参考方案（灰色）之间存在明显的增量：
这种行为将取决于帧速率，这可能会造成意想不到的后果。 我见过同事在一般的PC或笔记本电脑上运行未经优化的构建，并以中低帧频调整游戏，却发现在目标帧率下的感觉完全不同。 帧速率依赖性的典型症状是，在低帧速率下动态感觉松散和游动。 为了在低帧率和高帧率下都能保持一致，我们可以使弹簧以恒定的固定dt运行。 实际上，以小于最大允许时间步长的dt运行就足够了。
请注意，如果弹簧在物理更新中运行，则可能已经在使用固定的dts。 如果在这种情况下弹簧动力学存在问题，则可以以更高的频率（例如120Hz）运行弹簧动力学。
本节将介绍分步代码，以控制dt /仿真频率。
1.强制每个子步骤最大dt，允许混合dts
这将确保更新dt小于或等于规定的最大值。 它将以最大dt的步长前进，然后在最后一步中处理所有余数。
class Spring
{
// ... code from above ...
void UpdateWithSubsteps( const float dt )
{
const int maxUpdates = 8;
const float maxDt = 1.0f/60.0f;
int c = 0;
float timeLeft = dt;
while( timeLeft > 0.0f && c++ < maxUpdates)
{
float substepDt = min( timeLeft, maxDt );
UpdateSpring( substepDt );
timeLeft -= substepDt;
}
}
};
maxUpdates是可选的安全防护。 我在下面的注释部分中对此进行了一些考虑。
2.在每个子步骤中强制执行最大dt，保持统一的dts
该变量旨在计算小于最大值的统一更新dt。 使用此计算的更新dt模拟每个步骤。
class Spring
{
// ... code from above ...
void UpdateWithSubsteps( const float dt )
{
const float maxDt = 1.0f/60.0f;
float stepCountF = ceil(dt / maxDt);
float substepDt = dt / stepCountF;
int stepCount = (int)stepCountF;
const int maxStepCount = 8;
stepCount = min( stepCount, maxStepCount );
for( int i = 0; i < stepCount; i++ )
{
UpdateSpring( substepDt );
}
}
};
3.强制执行固定dt
这可能是这里概述的三个选项中最强大的，因为dt始终是固定的，并且更新逻辑将是确定性的和一致的（忽略影响浮点确定性的外部因素）。 但是，它也是最复杂的实现。
class Spring
{
// ... code from above ...
// the time the spring dynamics needs to be simulated forward
float m_timeToSimulate = 0.0f;
// the previous dynamic state
Vector3 m_prevPos, m_prevVel;
// the result state of the spring - the interpolated dynamic state that will be correct for the current frame
Vector3 m_currentFramePos, m_currentFrameVel;
void UpdateWithSubsteps( const float dt )
{
// we want to simulate forward by dt
m_timeToSimulate += dt;
const int maxUpdates = 8;
const float fixedDt = 1.0f/60.0f;
int c = 0;
while( m_timeToSimulate > 0.0f && c++ < maxUpdates )
{
// save pre-update dynamic state
m_prevPos = m_pos;
m_prevVel = m_vel;
// compute a new latest dynamic state. always the same dt.
UpdateSpring( fixedDt );
m_timeToSimulate -= fixedDt;
}
// at this point the spring has been updated beyond the target time. if we did not hit the max update limit,
// then m_timeToSimulate will be in (-fixedDt, 0]. Convert this to a lerp alpha:
float lerpAlpha = 1.0f + m_timeToSimulate / fixedDt;
// if we maxed out on updates above, we won't have managed to simulate up to the current time, in which
// case just use latest available state by clamping lerpAlpha
lerpAlpha = min( lerpAlpha, 1.0f );
// this is the public facing state of the spring that should be queried/used this frame
m_currentFramePos = lerp( m_prevPos, m_pos, lerpAlpha );
m_currentFrameVel = lerp( m_prevVel, m_vel, lerpAlpha );
}
// m_currentFramePos is the position that the spring should return to outside code
const Vector3& GetPos()
{
return m_currentFramePos;
}
// make sure we reset all of the state when poking data in
void SetPos( const Vector3& newPos )
{
m_currentFramePos = m_prevPos = m_pos = newPos;
}
};
杂项说明
maxUpdates是可选的安全防护。总的来说，这是一件好事–如果出现较长的帧故障，则下一个更新可能会进行许多spring更新以赶上实时传递的时间，这可能会进一步加剧性能问题（尽管不太可能在代码中添加简单的spring更新代码）以上会造成很大的损害）。
达到此更新限制后，更新将过早中止，并且不会一直模拟弹簧。看起来好像弹簧在缓慢运动。它也可能明显地抖动，因此会出现较大的帧故障。理想情况下，可以使所有框架都在预算资源之内。
必须小心确保将弹簧精确地更新为正确的量-dt向前模拟了弹簧的最终状态，否则系统将容易出现可见的抖动。确实很容易出错，我过去花了很多时间对其进行调试。我在这里详细介绍了抖动问题。
总结
这篇文章提出了一些对简单的弹簧实现的改进。 结果给出了一致的行为，在低帧率下具有顽健性，并且参数具有直观的含义并且易于调整和维护。
我建议您在不受物理更新影响的所有弹簧更新中添加子步骤，以确保弹簧行为在所有帧频下均有效。 由于它的简单性，第二个子步骤变体（均匀dts）是我的最爱。
References
[1] Tan J., Liu K., Turk G., Stable Proportional-Derivative Controllers,http://www.jie-tan.net/project/spd.pdf
[2] Bowles H., Implicit Spring, ShaderToy,https://www.shadertoy.com/view/MlB3Dm
Appendix A – Closed form solution
有一个封闭形式的弹簧动力学解析解，此处以灰色绘制：
此解决方案已在Wikipedia的the Universal oscillator equation section部分中记录。 在上面绘制的示例中，目标位置和速度不变，并且弹簧参数（UAF，DR）恒定，并且没有外力作用于弹簧。 在这种情况下，可以基于阻尼比选择正确的瞬态解，并使用起始条件（初始位置，速度）求解常数。 这是在ShaderToy（链接）中完成的。 在这种简单的场景中，不需要仿真-状态可以在将来的任意时间直接计算，而不必担心帧速率/ dts。
不幸的是，在实践中，我们很少遇到这样简单的场景，对于动力学使用封闭形式的表达式通常是不切实际的或无益的
Appendix B – Rigidbody control
如上所述，阻尼弹簧等效于PD控制，是控制物理刚体的有用工具。 下面的弹簧变量引用了刚体，并施加了将刚体移动到目标位置/速度的力：
class Spring
{
// target for spring, set as desired by user
Vector3 m_targetPos, m_targetVel;
// rigidbody
Rigidbody m_rb;
// spring param - Undamped Angular Frequency (UAF) - oscillation frequency
float m_UAF;
// spring param - Damping Ratio (DR)
float m_DR;
void PhysicsUpdate( const float dt )
{
// compute spring and damping forces
Vector3 posNextState = m_rb.Pos() + dt * m_rb.Vel();
Vector3 Fspring = m_UAF * m_UAF * (m_targetPos - posNextState);
Vector3 Fdamp = 2.0f * m_DR * m_UAF * (m_targetVel - m_rb.Vel());
// add acceleration
m_rb.Vel() += (Fspring + Fdamp) * dt / (1.0f + dt * 2.0f * m_DR * m_UAF);
}
};
现在，动态状态存在于刚体内。 不再需要用于位置和速度的成员变量，并且已将其删除。
第21行会显得很奇怪-我将加速度作为一种脉冲/速度变化来应用，而不是使用力。 使用力的结果将取决于RB的质量，这对我之前的项目没有兴趣。 部分原因是质量的影响并非微不足道–旋转惯量通常在轴上不对称，这是我不希望担心的。 这也部分是因为在先前的项目中，我们最终将质量设置为非物理值，以使物理“感觉”正确并保持稳定性。 因此，我发现将弹簧运动与质量隔离起来很有用。 这样，内容作者可以根据自己的喜好调整运动，而不必担心质量发生变化以后的运动变化。
缺点是必须调整UAF，以使其具有重量感。 如果不希望这样的话，则可以更改上述实现，以将Fspring和Fdamp作为力应用到刚体。
Acknowledgements
The work presented here is the outcome of collaboration between myself and my colleagues at Studio Gobo including Daniel Zimmermann, Chino Noris and Tom Williams.
**粗体** _斜体_ [链接](http://example.com) `代码` - 列表 > 引用。你还可以使用@来通知其他用户。