See also http://code.google.com/p/bullet/issues/detail?id=390 Added Demos/DX11ClothDemo (an OpenCL cloth demo will follow soon)
91 lines
2.8 KiB
Common Lisp
91 lines
2.8 KiB
Common Lisp
MSTRINGIFY(
|
|
|
|
/*#define float3 float4
|
|
float dot3(float3 a, float3 b)
|
|
{
|
|
return a.x*b.x + a.y*b.y + a.z*b.z;
|
|
}*/
|
|
|
|
float3 projectOnAxis( float3 v, float3 a )
|
|
{
|
|
return (a*dot(v, a));
|
|
}
|
|
|
|
__kernel void
|
|
ApplyForcesKernel(
|
|
const uint numNodes,
|
|
const float solverdt,
|
|
const float epsilon,
|
|
__global int * g_vertexClothIdentifier,
|
|
__global float4 * g_vertexNormal,
|
|
__global float * g_vertexArea,
|
|
__global float * g_vertexInverseMass,
|
|
__global float * g_clothLiftFactor,
|
|
__global float * g_clothDragFactor,
|
|
__global float4 * g_clothWindVelocity,
|
|
__global float4 * g_clothAcceleration,
|
|
__global float * g_clothMediumDensity,
|
|
__global float4 * g_vertexForceAccumulator,
|
|
__global float4 * g_vertexVelocity)
|
|
{
|
|
unsigned int nodeID = get_global_id(0);
|
|
if( nodeID < numNodes )
|
|
{
|
|
int clothId = g_vertexClothIdentifier[nodeID];
|
|
float nodeIM = g_vertexInverseMass[nodeID];
|
|
|
|
if( nodeIM > 0.0f )
|
|
{
|
|
float3 nodeV = g_vertexVelocity[nodeID].xyz;
|
|
float3 normal = g_vertexNormal[nodeID].xyz;
|
|
float area = g_vertexArea[nodeID];
|
|
float3 nodeF = g_vertexForceAccumulator[nodeID].xyz;
|
|
|
|
// Read per-cloth values
|
|
float3 clothAcceleration = g_clothAcceleration[clothId].xyz;
|
|
float3 clothWindVelocity = g_clothWindVelocity[clothId].xyz;
|
|
float liftFactor = g_clothLiftFactor[clothId];
|
|
float dragFactor = g_clothDragFactor[clothId];
|
|
float mediumDensity = g_clothMediumDensity[clothId];
|
|
|
|
// Apply the acceleration to the cloth rather than do this via a force
|
|
nodeV += (clothAcceleration*solverdt);
|
|
|
|
g_vertexVelocity[nodeID] = (float4)(nodeV, 0.f);
|
|
|
|
float3 relativeWindVelocity = nodeV - clothWindVelocity;
|
|
float relativeSpeedSquared = dot(relativeWindVelocity, relativeWindVelocity);
|
|
|
|
if( relativeSpeedSquared > epsilon )
|
|
{
|
|
// Correct direction of normal relative to wind direction and get dot product
|
|
normal = normal * (dot(normal, relativeWindVelocity) < 0 ? -1.f : 1.f);
|
|
float dvNormal = dot(normal, relativeWindVelocity);
|
|
if( dvNormal > 0 )
|
|
{
|
|
float3 force = (float3)(0.f, 0.f, 0.f);
|
|
float c0 = area * dvNormal * relativeSpeedSquared / 2.f;
|
|
float c1 = c0 * mediumDensity;
|
|
force += normal * (-c1 * liftFactor);
|
|
force += normalize(relativeWindVelocity)*(-c1 * dragFactor);
|
|
|
|
float dtim = solverdt * nodeIM;
|
|
float3 forceDTIM = force * dtim;
|
|
|
|
float3 nodeFPlusForce = nodeF + force;
|
|
|
|
// m_nodesf[i] -= ProjectOnAxis(m_nodesv[i], force.normalized())/dtim;
|
|
float3 nodeFMinus = nodeF - (projectOnAxis(nodeV, normalize(force))/dtim);
|
|
|
|
nodeF = nodeFPlusForce;
|
|
if( dot(forceDTIM, forceDTIM) > dot(nodeV, nodeV) )
|
|
nodeF = nodeFMinus;
|
|
|
|
g_vertexForceAccumulator[nodeID] = (float4)(nodeF, 0.0f);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
); |