在虚幻引擎中使用 Niagara 的平滑粒子流体动力学

如何解决在虚幻引擎中使用 Niagara 的平滑粒子流体动力学

这将是一个很长的帖子，抱歉，但我认为这是值得的，因为它非常复杂，我想很多其他人真的很想能够达到这种效果。这里还有一些关于 SPH 的其他问题，但没有一个与 Niagara 实施有关。我也在 Unreal Engine Answers 上发布了这个问题。

我一直在尝试在 Niagara 中复制流体模拟，如 Asher Zhu 所示：The Art of Illusion - Niagara Simulation Framework Overview。跳到 20:25 以获得我想要的效果。

看到他从渲染它的一些位（就我尚未得到的部分）中根本没有展示 Niagara 系统的任何部分，我遵循了此处的文章：link。

现在，我让它看起来或多或少像一种液体。然而，它看起来并不像亚瑟的。它相当不稳定，在爆炸然后安定下来之前，往往会在密度较高的区域停留几秒钟。它也永远不会发展任何深度。所有粒子，除非它们不规则地飞来飞去，否则都坐在地板上。另一个问题是碰撞——我看不出 Asher 是如何设法与环境发生如此干净的碰撞的。我的带符号距离场又大又圆且不均匀，粒子永远不会靠近墙壁。

下面的第四张图显示了它在到达第三张图后就爆炸了，第五张图是它最终稳定下来后的样子（以及粒子最终离墙壁多远）。最后一张图片显示它是完全平坦的（这不是盒子体积的问题；我已经测试过了）。

在这里很难显示 Niagara 系统中的所有内容，但关键是 HLSL 代码：

OutVelocity = Velocity;
OutPosition = Position;
Density = 0;
float Pressure = 0;

float smoothingRadius = 1.0f;
float restDensity = 0.2f;
float viscosity = 0.018f;
float gas = 500.0f;

const float3 gravity = float3(0,-98);
float pi = 3.141593;

int numParticles;
DirectReads.GetNumParticles(numParticles);

const float Poly6_constant = (315 / (64 * pi * pow(smoothingRadius,9)));
const float Spiky_constant = (-45 / (pi * pow(smoothingRadius,6)));

float3 forcePressure = float3(0,0);
float3 forceViscosity = float3(0,0);


#if GPU_SIMULATION

//Calculate the density of this particle based on the proximity of the other particles.
for (int i = 0; i < numParticles; ++i)
{
    bool myBool; //Temporary bool used to catch valid/invalid results for direct reads.

    float OtherMass;
    DirectReads.GetFloatByIndex<Attribute="Mass">(i,myBool,OtherMass);
    float3 OtherPosition;
    DirectReads.GetVectorByIndex<Attribute="Position">(i,OtherPosition);

    // Calculate the distance and direction between the target Particle and itself
    float distanceBetween = distance(OtherPosition,OutPosition);
    
    if (distanceBetween < smoothingRadius)
    {
        Density += OtherMass * Poly6_constant * pow(smoothingRadius - distanceBetween,3);
    }
}

//Avoid negative pressure by clamping density to reference value
Density = max(restDensity,Density);

//Calculate pressure
Pressure = gas * (Density - restDensity);

//Calculate the forces.
for (int i = 0; i < numParticles; ++i)
{
    if (i != InstanceId) //Only calculate the pressure-based force and Laplacian smoothing function if the other particle is not the current particle.)
    {
        bool myBool; //Temporary bool used to catch valid/invalid results for direct reads.

        float OtherMass;
        DirectReads.GetFloatByIndex<Attribute="Mass">(i,OtherMass);
        float OtherDensity;
        DirectReads.GetFloatByIndex<Attribute="Density">(i,OtherDensity);
        float3 OtherPosition;
        DirectReads.GetVectorByIndex<Attribute="Position">(i,OtherPosition);
        float3 OtherVelocity;
        DirectReads.GetVectorByIndex<Attribute="Velocity">(i,OtherVelocity);

        float3 direction = OutPosition - OtherPosition;
        float3 normalisedVector = normalize(direction);
        float distanceBetween = distance(OtherPosition,OutPosition);

        if (distanceBetween > 0 && distanceBetween < smoothingRadius) //distanceBetween must be >0 to avoide a div0 error.
        {
            float OtherPressure = gas * (OtherDensity - restDensity);

            //Calculate particle pressure.
            forcePressure += -1 * Mass * normalisedVector * (Pressure + OtherPressure) / (2 * Density * OtherDensity) * Spiky_constant * pow(smoothingRadius - distanceBetween,2);

            //Viscosity-based force computation with Laplacian smoothing function (W).
               const float W = -(pow(distanceBetween,3) / (2 * pow(smoothingRadius,3))) + (pow(distanceBetween,2) / pow(smoothingRadius,2)) + (smoothingRadius / (2 * distanceBetween)) - 1;
            forceViscosity += viscosity * (OtherMass / Mass) * (1 / OtherDensity) * (OtherVelocity - Velocity) * W * normalisedVector;
            //forceViscosity += viscosity * (OtherMass / Mass) * (1 / OtherDensity) * (OtherVelocity - Velocity) * (45 / (pi * pow(smoothingRadius,6))) * (smoothingRadius - distanceBetween);
        }
    }
}

OutVelocity += DeltaTime * ((forcePressure + forceViscosity) / Density);
OutPosition += DeltaTime * OutVelocity;
#endif

此代码对系统中的所有其他粒子执行两次循环，一次计算压力，另一次计算力。然后它输出速度和位置。就像我上面链接的文章一样，就像我看到的其他一些东西一样。然而，它的行为根本不像那些资源中显示的那样。

我没有应用任何基于网格的优化。为此，我将仅应用 UE 的内容示例项目中的 PBD 示例中使用的网格优化。但就目前而言，这是一个并不真正需要的额外复杂性。即使没有它，它也能在数千个粒子下运行良好。

我查看了一些资源（文章、视频和学术研究论文），并花了两周时间进行试验，包括对代码顶部的值进行反复试验。我显然错过了一些关键的东西。它可以是什么？我现在很沮丧，任何帮助都将不胜感激。