基于距离约束

每条边对应一个距离约束

float compress = 1.0f, stretch = 1.0f;
for(uint32_t i = 0; i < cloak.indices.size(); i = i + 3)
{
	size_t indexA = cloak.indices[i];
	size_t indexB = cloak.indices[i+1];
	size_t indexC = cloak.indices[i+2];

	PiratePhysics::Rigidbody &a = *(pbd.rigidbodies[2+indexA]);
	PiratePhysics::Rigidbody &b = *(pbd.rigidbodies[2+indexB]);
	PiratePhysics::Rigidbody &c = *(pbd.rigidbodies[2+indexC]);

	float distAB = (a.position - b.position).norm();
	pbd.constraints.emplace_back(std::make_shared<Stretching>(a, b, distAB, compress, stretch));

	float distAC = (a.position - c.position).norm();
	pbd.constraints.emplace_back(std::make_shared<Stretching>(a, c, distAC, compress, stretch));

	float distBC = (b.position - c.position).norm();
	pbd.constraints.emplace_back(std::make_shared<Stretching>(b, c, distBC, compress, stretch));
}

空气动力学模型

$\mathbf{u}$ 风力速度场， $\mathbf{F}_{L}$ 升力， $\mathbf{F}_D$ ， $\mathbf{\hat{n}}$ 法线， $\mathbf{V}_{rel}$ 粒子相对速度， $q$ 相对速度方向与法线方向的夹角
布料的每个三角面的三个粒子 $a,b,c$ ，则：
$\mathbf{V}_{rel} = \frac{1}{3}(\mathbf{v}_a+\mathbf{v}_b+\mathbf{v}_c) - \mathbf{u}$
法线方向的选取使得：
$\mathbf{V}_{rel} \cdot \mathbf{\hat{n}} > 0$
升力：
$\mathbf{F}_{L} = \frac{1}{2}C_L\rho|\mathbf{V}_{rel}|^2A\cos\theta \mathbf{\hat{n}}\times \mathbf{V}_{rel} \times \mathbf{V}_{rel}$
阻力：
$\mathbf{F}_{D} = \frac{1}{2}C_D\rho|\mathbf{V}_{rel}|^2A (\mathbf{\hat{n}}\cdot \mathbf{V}_{rel}) (-\mathbf{V}_{rel})$
$C_L$ 升力系数， $C_D$ 阻力系数， $\rho$ 空气密度

Vector3f u{2.f, 0.f, 0.f};
float drag = 1.f;
float lift = 1.f;

for(size_t i = 0; i < cloak.indices.size(); i = i + 3)
{
	size_t indexA = cloak.indices[i];
	size_t indexB = cloak.indices[i+1];
	size_t indexC = cloak.indices[i+2];
	PiratePhysics::Particle &a = *dynamic_cast<Particle*>(pbd.rigidbodies[2+indexA].get());
	PiratePhysics::Particle &b = *dynamic_cast<Particle*>(pbd.rigidbodies[2+indexB].get());
	PiratePhysics::Particle &c = *dynamic_cast<Particle*>(pbd.rigidbodies[2+indexC].get());

	Vector3f AB = a.positionPrime - b.positionPrime;
	Vector3f AC = a.positionPrime - c.positionPrime;
	Vector3f ABC = AB.cross(AC);
	float area = ABC.norm() / 2;

	Vector3f vRel = (a.velocity + b.velocity + c.velocity)/3 - u;
	Vector3f n = ABC.normalized();
	n = vRel.dot(n) > 0 ? n : -n;

	Vector3f fd = drag * vRel.squaredNorm() * area * vRel.dot(n) * -vRel;
	float cosq = vRel.normalized().cross(n).norm();
	Vector3f fl = lift * vRel.squaredNorm() * area * cosq * (n.cross(vRel).cross(vRel));

	a.velocity += dt * a.invMass *(fd+fl);
	b.velocity += dt * b.invMass *(fd+fl);
	c.velocity += dt * c.invMass *(fd+fl);
}

Reference

Keckeisen, Michael, et al. "Modelling effects of wind fields in cloth animations." (2004).

布料模拟

基于距离约束

空气动力学模型

Reference