转自
前言
声明,此 NeHe OpenGL教程系列文章由51博客翻译(2010-08-19),本博客为转载并稍加整理与修改。对NeHe的OpenGL管线教程的编写,以及yarn的翻译整理表示感谢。
NeHe OpenGL第四十课:绳子的模拟
绳子的模拟:
怎样模拟一根绳子呢,把它想象成一个个紧密排列的点,怎么样有了思路了吧,在这一课你你将学会怎样建模,简单吧,你能模拟更多。 绳索模拟
在这个教程里我们将模拟一段绳索,我们是在39课的基础上进行的。
在物理模拟中,我们必须设置各个物理量,就像它们在自然界中的行为一样。模拟中的运动并不一定和自然界相同,我们使用的运动模型,必须和我们需要模拟的目的有关,目的决定了它的精确度。要知道我们的目标不
是模拟原子和分子,也不是模拟成千上万的粒子系。首先我们需要确定我们模拟的目标,才能创建我们的物理模型。它和下面内容相关:
1. 运动的数学表示 2. 执行模拟的计算机的速度
1. 运动的数学表示:
这个问题决定了我们使用何种数学方程来模拟运动,使用经典力学还是量子力学。
2. 执行模拟的计算机的速度:
计算机的速度决定了我们可以模拟的精度。
设计绳索的物理模型:
我们在经典力学和高于500Mhz的计算机上模拟这个问题。首先我们需要设定需要的精度,我们使用一系列互相用弹簧连接的质点来模拟绳索,精度决定了我们用多少个点来模拟,当然越多越精确。在下面我决定用50或
100个点来模拟绳子一段3或4m长的绳子,换句话说,我们的模拟精度就是3到8厘米。
设计运动模型:
在绳子中,施加给各个质点的力来自于自身的质量和相连的内力(参见大学里的普通力学)。如下我们用"O"表示质点,“—”表示连接质点的弹簧。 O----O----O----O 1 2 3 4
弹簧的力学公式如下:
力 = -k * x k: 弹性系数 x: 相距平衡位置的位移
上面的公式说明,如果相邻点的距离为平衡距离,那么它们不受到任何力的作用。如果我们设置平衡位置为5cm,那么100个点的绳子长5m。如果相连质点之间的位置小于5cm,它们受到排斥力。
上面的公式只是一个基础,现在我们可以加上摩擦力,如果没有这项,那么绳子将永远动下去。
弹簧类:
这个类包含相连接的两个物体,它们之间具有作用力。 class Spring { public: Mass* mass1; // 质点1 Mass* mass2; // 质点2
float springConstant; // 弹性系数 float springLength; //弹簧长度 float frictionConstant; //摩擦系数
Spring(Mass* mass1, Mass* mass2, // 构造函数 float springConstant, float springLength, float frictionConstant) { this->springConstant = springConstant; this->springLength = springLength; this->frictionConstant = frictionConstant;
this->mass1 = mass1; this->mass2 = mass2; }
void solve() // 计算各个物体的受力 { Vector3D springVector = mass1->pos - mass2->pos; float r = springVector.length(); // 计算两个物体之间的距离
Vector3D force; if (r != 0) // 计算力 force += -(springVector / r) * (r - springLength) * springConstant; ...
force += -(mass1->vel - mass2->vel) * frictionConstant; // 加上摩擦力 mass1->applyForce(force); // 给物体1施加力 mass2->applyForce(-force); // 给物体2施加力 } 下面我们把绳子钉在墙上,所以我们的模拟就多了一个万有引力,空气摩擦力。万有引力的公式如下:
力 = (重力加速度) * 质量
万有引力会作用在每一个质点上,地面也会给每个物体一个作用力。在我们的模型中将考虑绳子和地面之间的接触,地面给绳子向上的力,并提供摩擦力。
设置模拟的初始值
现在我们已经设置好模拟环境了,长度单位是m,时间单位是秒,质量单位是kg。
为了设置初始值,我们必须提供供模拟开始的参数。我们定义一下参数:
1. 重力加速度: 9.81 m/s/s 垂直向下 2. 质点个数: 80 3. 相连质点的距离: 5 cm (0.05 meters) 4. 质量: 50 克(0.05 kg) 5. 绳子开始处于垂直状态
下面计算绳子受到的力
f = (绳子质量) * (重力加速度) = (4 kg) * (9.81) ~= 40 N
弹簧必须平衡这个力 40 N,它伸长1cm,计算弹性系数:
合力= -k * x = -k * 0.01 m
合力应该为0 :
40 N + (-k * 0.01 meters) = 0
弹性系数 k 为:
k = 4000 N / m
设置弹簧的摩擦系数:
springFrictionConstant = 0.2 N/(m/s)
下面我们看看这个绳索类:
1. virtual void init() ---> 重置力
2. virtual void solve() ---> 计算各个质点的力
3. virtual void simulate(float dt) ---> 模拟一次
4. virtual void operate(float dt) ---> 执行一次操作
绳索类如下所示 : class RopeSimulation : public Simulation //绳索类 { public: Spring** springs; // 弹簧类结构的数组的指针
Vector3D gravitation; // 万有引力
Vector3D ropeConnectionPos; // 绳索的连接点 Vector3D ropeConnectionVel; //连接点的速度,我们使用这个移动绳子
float groundRepulsionConstant; //地面的反作用力 float groundFrictionConstant; //地面的摩擦系数 float groundAbsorptionConstant; //地面的缓冲力 float groundHeight; //地面高度
float airFrictionConstant; //空气的摩擦系数 下面是它的构造函数
RopeSimulation( int numOfMasses, float m, float springConstant, float springLength, float springFrictionConstant, Vector3D gravitation, float airFrictionConstant, float groundRepulsionConstant, float groundFrictionConstant, float groundAbsorptionConstant, float groundHeight ) : Simulation(numOfMasses, m) { this->gravitation = gravitation; this->airFrictionConstant = airFrictionConstant;
this->groundFrictionConstant = groundFrictionConstant; this->groundRepulsionConstant = groundRepulsionConstant; this->groundAbsorptionConstant = groundAbsorptionConstant; this->groundHeight = groundHeight;
for (int a = 0; a < numOfMasses; ++a) // 设置质点位置 { masses[a]->pos.x = a * springLength; masses[a]->pos.y = 0; masses[a]->pos.z = 0; }
springs = new Spring*[numOfMasses - 1]; for (a = 0; a < numOfMasses - 1; ++a) //创建各个质点之间的模拟弹簧 { springs[a] = new Spring(masses[a], masses[a + 1], springConstant, springLength, springFrictionConstant); } }
计算施加给各个质点的力
void solve() // 计算施加给各个质点的力 { for (int a = 0; a < numOfMasses - 1; ++a) // 弹簧施加给各个物体的力 { springs[a]->solve(); }
for (a = 0; a < numOfMasses; ++a) // 计算各个物体受到的其它的力 { masses[a]->applyForce(gravitation * masses[a]->m); // 万有引力 // 空气的摩擦力 masses[a]->applyForce(-masses[a]->vel * airFrictionConstant);
if (masses[a]->pos.y < groundHeight) // 计算地面对质点的作用 { Vector3D v;
v = masses[a]->vel; // 返回速度 v.y = 0; // y方向的速度为0
// 计算地面给质点的力 masses[a]->applyForce(-v * groundFrictionConstant);
v = masses[a]->vel; v.x = 0; v.z = 0;
if (v.y < 0) // 计算地面的缓冲力
masses[a]->applyForce(-v * groundAbsorptionConstant); // 计算地面的反作用力 Vector3D force = Vector3D(0, groundRepulsionConstant, 0) * (groundHeight - masses[a]->pos.y);
masses[a]->applyForce(force); // 施加地面对质点的力 } } }
下面的代码完成整个模拟过程
void simulate(float dt) // 模拟一次 { Simulation::simulate(dt); // 调用基类的模拟函数
ropeConnectionPos += ropeConnectionVel * dt; // 计算绳子的连接点
if (ropeConnectionPos.y < groundHeight) { ropeConnectionPos.y = groundHeight; ropeConnectionVel.y = 0; }
masses[0]->pos = ropeConnectionPos; // 更新绳子的连接点和速度 masses[0]->vel = ropeConnectionVel; }
void setRopeConnectionVel(Vector3D ropeConnectionVel) { this->ropeConnectionVel = ropeConnectionVel; }
有了上面的类,我们可以很方便的模拟绳子,代码如下: RopeSimulation* ropeSimulation = new RopeSimulation( 80, // 80 质点 0.05f, // 每个质点50g 10000.0f, // 弹性系数 0.05f, // 质点之间的距离 0.2f, // 弹簧的内摩擦力 Vector3D(0, -9.81f, 0), // 万有引力 0.02f, // 空气摩擦力 100.0f, // 地面反作用系数 0.2f, // 地面摩擦系数 2.0f, // 地面缓冲系数 -1.5f); // 地面高度
下面的代码在程序中执行绳子的模拟 float dt = milliseconds / 1000.0f; // 经过的秒数
float maxPossible_dt = 0.002f; // 模拟间隔
int numOfIterations = (int)(dt / maxPossible_dt) + 1; // 模拟次数 if (numOfIterations != 0) dt = dt / numOfIterations;
for (int a = 0; a < numOfIterations; ++a) // 执行模拟 ropeSimulation->operate(dt); 原文及其个版本源代码下载: