去掉obi,使用自写绳索

This commit is contained in:
2026-02-23 20:51:03 +08:00
parent cb636f862d
commit 91e2309eeb
2011 changed files with 2593 additions and 190578 deletions

View File

@@ -1,553 +0,0 @@
#if (OBI_BURST && OBI_MATHEMATICS && OBI_COLLECTIONS)
using Unity.Burst;
using Unity.Collections;
using Unity.Collections.LowLevel.Unsafe;
using Unity.Jobs;
using Unity.Mathematics;
namespace Obi
{
public class BurstDensityConstraints : BurstConstraintsImpl<BurstDensityConstraintsBatch>
{
public NativeList<int> fluidParticles;
public BurstDensityConstraints(BurstSolverImpl solver) : base(solver, Oni.ConstraintType.Density)
{
fluidParticles = new NativeList<int>(Allocator.Persistent);
}
public override IConstraintsBatchImpl CreateConstraintsBatch()
{
var dataBatch = new BurstDensityConstraintsBatch(this);
batches.Add(dataBatch);
return dataBatch;
}
public override void Dispose()
{
fluidParticles.Dispose();
}
public override void RemoveBatch(IConstraintsBatchImpl batch)
{
batches.Remove(batch as BurstDensityConstraintsBatch);
batch.Destroy();
}
protected override JobHandle EvaluateSequential(JobHandle inputDeps, float stepTime, float substepTime, int steps, float timeLeft)
{
return EvaluateParallel(inputDeps, stepTime, substepTime, steps, timeLeft);
}
protected override JobHandle EvaluateParallel(JobHandle inputDeps, float stepTime, float substepTime, int steps, float timeLeft)
{
inputDeps = UpdateInteractions(inputDeps);
// evaluate all batches as a chain of dependencies:
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].Evaluate(inputDeps, stepTime, substepTime, steps, timeLeft);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
// calculate per-particle density lambdas:
inputDeps = CalculateLambdas(inputDeps, substepTime);
// calculate viscosity/vorticity:
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].ViscosityAndVorticity(inputDeps);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
// apply viscosity/vorticity positional deltas:
var app = new ApplyPositionDeltasJob()
{
fluidParticles = fluidParticles,
positions = m_Solver.positions,
deltas = m_Solver.positionDeltas,
counts = m_Solver.positionConstraintCounts,
anisotropies = m_Solver.anisotropies,
normals = m_Solver.normals,
fluidData = m_Solver.fluidData,
matchingRotations = m_Solver.orientationDeltas,
linearFromAngular = m_Solver.restPositions,
};
inputDeps = app.Schedule(fluidParticles.Length, 64, inputDeps);
// apply density positional deltas:
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].Apply(inputDeps, substepTime);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
return inputDeps;
}
public JobHandle CalculateVelocityCorrections(JobHandle inputDeps, float deltaTime)
{
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].CalculateNormals(inputDeps, deltaTime);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
return inputDeps;
}
public JobHandle ApplyVelocityCorrections(JobHandle inputDeps, float deltaTime)
{
inputDeps = ApplyAtmosphere(inputDeps, deltaTime);
m_Solver.ScheduleBatchedJobsIfNeeded();
return inputDeps;
}
public JobHandle CalculateAnisotropyLaplacianSmoothing(JobHandle inputDeps)
{
// if the constraints are deactivated or we need no anisotropy:
if (((BurstSolverImpl)solver).abstraction.parameters.maxAnisotropy <= 1)
return inputDeps;
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].AccumulateSmoothPositions(inputDeps);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
inputDeps = AverageSmoothPositions(inputDeps);
for (int i = 0; i < batches.Count; ++i)
{
if (batches[i].enabled)
{
inputDeps = batches[i].AccumulateAnisotropy(inputDeps);
m_Solver.ScheduleBatchedJobsIfNeeded();
}
}
return AverageAnisotropy(inputDeps);
}
private JobHandle UpdateInteractions(JobHandle inputDeps)
{
// clear existing fluid data:
var clearData = new ClearFluidDataJob()
{
fluidParticles = fluidParticles,
fluidData = m_Solver.fluidData,
massCenters = m_Solver.normals,
prevMassCenters = m_Solver.renderablePositions,
moments = m_Solver.anisotropies
};
inputDeps = clearData.Schedule(fluidParticles.Length, 64, inputDeps);
// update fluid interactions:
var updateInteractions = new UpdateInteractionsJob()
{
pairs = m_Solver.fluidInteractions,
positions = m_Solver.positions,
fluidMaterials = m_Solver.fluidMaterials,
densityKernel = new Poly6Kernel(((BurstSolverImpl)solver).abstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
gradientKernel = new SpikyKernel(((BurstSolverImpl)solver).abstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
};
return updateInteractions.Schedule(((BurstSolverImpl)solver).fluidInteractions.Length, 64, inputDeps);
}
private JobHandle CalculateLambdas(JobHandle inputDeps, float deltaTime)
{
// calculate lagrange multipliers:
var calculateLambdas = new CalculateLambdasJob
{
fluidParticles = fluidParticles,
positions = m_Solver.positions,
prevPositions = m_Solver.prevPositions,
matchingRotations = m_Solver.restPositions.Reinterpret<quaternion>(),
principalRadii = m_Solver.principalRadii,
fluidMaterials = m_Solver.fluidMaterials,
densityKernel = new Poly6Kernel(m_Solver.abstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
gradientKernel = new SpikyKernel(m_Solver.abstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
fluidData = m_Solver.fluidData,
massCenters = m_Solver.normals,
prevMassCenters = m_Solver.renderablePositions,
moments = m_Solver.anisotropies,
deltas = m_Solver.positionDeltas,
counts = m_Solver.positionConstraintCounts,
solverParams = m_Solver.abstraction.parameters
};
return calculateLambdas.Schedule(fluidParticles.Length,64,inputDeps);
}
private JobHandle ApplyAtmosphere(JobHandle inputDeps, float deltaTime)
{
var conf = new ApplyAtmosphereJob
{
fluidParticles = fluidParticles,
wind = m_Solver.wind,
fluidInterface = m_Solver.fluidInterface,
fluidMaterials2 = m_Solver.fluidMaterials2,
principalRadii = m_Solver.principalRadii,
normals = m_Solver.normals,
fluidData = m_Solver.fluidData,
velocities = m_Solver.velocities,
angularVelocities = m_Solver.angularVelocities,
vorticity = m_Solver.restOrientations.Reinterpret<float4>(),
vorticityAccelerations = m_Solver.orientationDeltas.Reinterpret<float4>(),
linearAccelerations = m_Solver.positionDeltas,
linearFromAngular = m_Solver.restPositions,
angularDiffusion = m_Solver.anisotropies,
positions = m_Solver.positions,
prevPositions = m_Solver.prevPositions,
dt = deltaTime,
solverParams = m_Solver.abstraction.parameters
};
return conf.Schedule(fluidParticles.Length, 64, inputDeps);
}
private JobHandle AverageSmoothPositions(JobHandle inputDeps)
{
var average = new AverageSmoothPositionsJob()
{
fluidParticles = fluidParticles,
renderablePositions = m_Solver.renderablePositions,
anisotropies = m_Solver.anisotropies
};
return average.Schedule(fluidParticles.Length, 64, inputDeps);
}
private JobHandle AverageAnisotropy(JobHandle inputDeps)
{
var average = new AverageAnisotropyJob()
{
fluidParticles = fluidParticles,
renderablePositions = m_Solver.renderablePositions,
renderableOrientations = m_Solver.renderableOrientations,
principalRadii = m_Solver.principalRadii,
anisotropies = m_Solver.anisotropies,
maxAnisotropy = m_Solver.abstraction.parameters.maxAnisotropy,
renderableRadii = m_Solver.renderableRadii,
fluidData = m_Solver.fluidData,
life = m_Solver.life,
solverParams = m_Solver.abstraction.parameters
};
return average.Schedule(fluidParticles.Length, 64, inputDeps);
}
[BurstCompile]
public struct ClearFluidDataJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> massCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> prevMassCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> moments;
public void Execute(int i)
{
int p = fluidParticles[i];
fluidData[p] = float4.zero;
massCenters[p] = float4.zero;
prevMassCenters[p] = float4.zero;
moments[p] = float4x4.zero;
}
}
[BurstCompile]
public struct UpdateInteractionsJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> positions;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public Poly6Kernel densityKernel;
[ReadOnly] public SpikyKernel gradientKernel;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<FluidInteraction> pairs;
[ReadOnly] public BatchData batchData;
public void Execute(int i)
{
var pair = pairs[i];
// calculate normalized gradient vector:
pair.gradient = new float4((positions[pair.particleA] - positions[pair.particleB]).xyz,0);
float distance = math.length(pair.gradient);
pair.gradient /= distance + math.FLT_MIN_NORMAL;
// calculate and store average density and gradient kernels:
pair.avgKernel = (densityKernel.W(distance, fluidMaterials[pair.particleA].x) +
densityKernel.W(distance, fluidMaterials[pair.particleB].x)) * 0.5f;
pair.avgGradient = (gradientKernel.W(distance, fluidMaterials[pair.particleA].x) +
gradientKernel.W(distance, fluidMaterials[pair.particleB].x)) * 0.5f;
pairs[i] = pair;
}
}
[BurstCompile]
public struct CalculateLambdasJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[ReadOnly] public NativeArray<float4> positions;
[ReadOnly] public NativeArray<float4> prevPositions;
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public Poly6Kernel densityKernel;
[ReadOnly] public SpikyKernel gradientKernel;
[ReadOnly] public Oni.SolverParameters solverParams;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> massCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> prevMassCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> moments;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<quaternion> matchingRotations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> deltas;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<int> counts;
public void Execute(int p)
{
int i = fluidParticles[p];
float restVolume = math.pow(principalRadii[i].x * 2, 3 - (int)solverParams.mode);
float4 data = fluidData[i];
float grad = restVolume * gradientKernel.W(0, fluidMaterials[i].x);
// self particle contribution to density, gradient and mass centers:
data += new float4(densityKernel.W(0, fluidMaterials[i].x), 0, grad, grad * grad + data[2] * data[2]);
massCenters[i] += new float4(positions[i].xyz, 1) / positions[i].w;
prevMassCenters[i] += new float4(prevPositions[i].xyz, 1) / positions[i].w;
// usually, we'd weight density by mass (density contrast formulation) by dividing by invMass. Then, multiply by invMass when
// calculating the state equation (density / restDensity - 1, restDensity = mass / volume, so density * invMass * restVolume - 1
// We end up with density / invMass * invMass * restVolume - 1, invMass cancels out.
float constraint = math.max(0, data[0] * restVolume - 1) * fluidMaterials[i].w;
// calculate lambda:
data[1] = -constraint / (positions[i].w * data[3] + math.FLT_MIN_NORMAL);
fluidData[i] = data;
// get total neighborhood mass:
float M = massCenters[i][3];
massCenters[i] /= massCenters[i][3];
prevMassCenters[i] /= prevMassCenters[i][3];
// update moments:
moments[i] += (BurstMath.multrnsp4(positions[i], prevPositions[i]) + float4x4.identity * math.pow(principalRadii[i].x, 2) * 0.001f) / positions[i].w;
moments[i] -= M * BurstMath.multrnsp4(massCenters[i], prevMassCenters[i]);
// extract neighborhood orientation delta:
matchingRotations[i] = BurstMath.ExtractRotation(moments[i], quaternion.identity, 5);
// viscosity and vorticity:
float4 viscGoal = new float4(massCenters[i].xyz + math.rotate(matchingRotations[i], (prevPositions[i] - prevMassCenters[i]).xyz), 0);
deltas[i] += (viscGoal - positions[i]) * fluidMaterials[i].z;
counts[i]++;
}
}
[BurstCompile]
public struct ApplyPositionDeltasJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> positions;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> deltas;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<int> counts;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> normals;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> anisotropies;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<quaternion> matchingRotations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> linearFromAngular;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
public void Execute(int p)
{
int i = fluidParticles[p];
if (counts[i] > 0)
{
positions[i] += new float4(deltas[i].xyz,0) / counts[i];
deltas[i] = float4.zero;
counts[i] = 0;
}
normals[i] = float4.zero;
anisotropies[i] = float4x4.zero;
linearFromAngular[i] = float4.zero;
matchingRotations[i] = new quaternion(0, 0, 0, 0);
// zero out fluidData.z in preparation to accumulate relative velocity.
float4 data = fluidData[i];
data.z = 0;
fluidData[i] = data;
}
}
[BurstCompile]
public struct ApplyAtmosphereJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[ReadOnly] public NativeArray<float4> wind;
[ReadOnly] public NativeArray<float4> fluidInterface;
[ReadOnly] public NativeArray<float4> fluidMaterials2;
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public NativeArray<float4> normals;
[ReadOnly] public NativeArray<float4> fluidData;
[ReadOnly] public NativeArray<float4> linearFromAngular;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> positions;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> prevPositions;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> linearAccelerations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> vorticityAccelerations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> vorticity;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> angularDiffusion;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> angularVelocities;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> velocities;
[ReadOnly] public float dt;
[ReadOnly] public Oni.SolverParameters solverParams;
public void Execute(int p)
{
int i = fluidParticles[p];
float restVolume = math.pow(principalRadii[i].x * 2, 3 - (int)solverParams.mode);
//atmospheric drag:
float4 velocityDiff = velocities[i] - wind[i];
// particles near the surface should experience drag:
velocities[i] -= fluidInterface[i].x * velocityDiff * math.max(0, 1 - fluidData[i][0] * restVolume) * dt;
// ambient pressure:
velocities[i] += fluidInterface[i].y * normals[i] * dt;
// angular accel due to baroclinity:
angularVelocities[i] += new float4(fluidMaterials2[i].z * math.cross(-normals[i].xyz, -velocityDiff.xyz), 0) * dt;
angularVelocities[i] -= fluidMaterials2[i].w * angularDiffusion[i].c0;
// micropolar vorticity:
velocities[i] += fluidMaterials2[i].x * linearAccelerations[i] * dt;
vorticity[i] += fluidMaterials2[i].x * (vorticityAccelerations[i] * 0.5f - vorticity[i]) * dt;
vorticity[i] -= fluidMaterials2[i].y * angularDiffusion[i].c1;
linearAccelerations[i] = float4.zero;
vorticityAccelerations[i] = float4.zero;
angularDiffusion[i] = float4x4.zero;
// we want to add together linear and angular velocity fields and use result to advect particles without modifying either field:
positions[i] += new float4(linearFromAngular[i].xyz * dt,0);
prevPositions[i] += new float4(linearFromAngular[i].xyz * dt, 0);
}
}
[BurstCompile]
public struct AverageSmoothPositionsJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[ReadOnly] public NativeArray<float4> renderablePositions;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4x4> anisotropies;
public void Execute(int p)
{
int i = fluidParticles[p];
var smoothPos = anisotropies[i];
if (smoothPos.c3.w > 0)
smoothPos.c3 /= smoothPos.c3.w;
else
smoothPos.c3.xyz = renderablePositions[i].xyz;
anisotropies[i] = smoothPos;
}
}
[BurstCompile]
public struct AverageAnisotropyJob : IJobParallelFor
{
[ReadOnly] public NativeList<int> fluidParticles;
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public float maxAnisotropy;
[ReadOnly] public NativeArray<float4x4> anisotropies;
[ReadOnly] public NativeArray<float> life;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> renderablePositions;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<quaternion> renderableOrientations;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> renderableRadii;
[ReadOnly] public Oni.SolverParameters solverParams;
public void Execute(int p)
{
int i = fluidParticles[p];
if (anisotropies[i].c3.w > 0 && (anisotropies[i].c0[0] + anisotropies[i].c1[1] + anisotropies[i].c2[2]) > 0.01f)
{
float3 singularValues;
float3x3 u;
BurstMath.EigenSolve(math.float3x3(anisotropies[i] / anisotropies[i].c3.w), out singularValues, out u);
float max = singularValues[0];
float3 s = math.max(singularValues,new float3(max / maxAnisotropy)) / max * principalRadii[i].x;
renderableOrientations[i] = quaternion.LookRotationSafe(u.c2,u.c1);
renderableRadii[i] = new float4(s.xyz,1);
}
else
{
float radius = principalRadii[i].x / maxAnisotropy;
renderableOrientations[i] = quaternion.identity;
renderableRadii[i] = new float4(radius,radius,radius,1);
float4 data = fluidData[i];
data.x = 1 / math.pow(math.abs(radius * 2), 3 - (int)solverParams.mode); // normal volume of an isolated particle.
fluidData[i] = data;
}
renderablePositions[i] = math.lerp(renderablePositions[i], anisotropies[i].c3, math.min((maxAnisotropy - 1)/3.0f,1));
// inactive particles have radii.w == 0, set it right away for particles killed during this frame
// to keep them from being rendered during this frame instead of waiting to do it at the start of next sim step:
float4 radii = renderableRadii[i];
radii.w = life[i] <= 0 ? 0 : radii.w;
renderableRadii[i] = radii;
}
}
}
}
#endif

View File

@@ -1,11 +0,0 @@
fileFormatVersion: 2
guid: df96ac4db14c846ea99a6ebc5771098f
MonoImporter:
externalObjects: {}
serializedVersion: 2
defaultReferences: []
executionOrder: 0
icon: {instanceID: 0}
userData:
assetBundleName:
assetBundleVariant:

View File

@@ -1,477 +0,0 @@
#if (OBI_BURST && OBI_MATHEMATICS && OBI_COLLECTIONS)
using Unity.Jobs;
using Unity.Collections;
using Unity.Collections.LowLevel.Unsafe;
using Unity.Mathematics;
using Unity.Burst;
namespace Obi
{
public class BurstDensityConstraintsBatch : BurstConstraintsBatchImpl, IDensityConstraintsBatchImpl
{
public BatchData batchData;
public BurstDensityConstraintsBatch(BurstDensityConstraints constraints)
{
m_Constraints = constraints;
m_ConstraintType = Oni.ConstraintType.Density;
}
public override JobHandle Initialize(JobHandle inputDeps, float stepTime, float substepTime, int steps, float timeLeft)
{
return inputDeps;
}
public override JobHandle Evaluate(JobHandle inputDeps, float stepTime, float substepTime, int steps, float timeLeft)
{
// update densities and gradients:
var updateDensities = new UpdateDensitiesJob()
{
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
positions = solverImplementation.positions,
prevPositions = solverImplementation.prevPositions,
principalRadii = solverImplementation.principalRadii,
fluidMaterials = solverImplementation.fluidMaterials,
fluidData = solverImplementation.fluidData,
moments = solverImplementation.anisotropies,
massCenters = solverImplementation.normals,
prevMassCenters = solverImplementation.renderablePositions,
densityKernel = new Poly6Kernel(solverAbstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
batchData = batchData,
solverParams = solverAbstraction.parameters
};
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
return updateDensities.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
public override JobHandle Apply(JobHandle inputDeps, float substepTime)
{
var parameters = solverAbstraction.GetConstraintParameters(m_ConstraintType);
// update densities and gradients:
var apply = new ApplyDensityConstraintsJob()
{
principalRadii = solverImplementation.principalRadii,
fluidMaterials = solverImplementation.fluidMaterials,
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
densityKernel = new Poly6Kernel(solverAbstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
positions = solverImplementation.positions,
fluidData = solverImplementation.fluidData,
batchData = batchData,
solverParams = solverAbstraction.parameters,
sorFactor = parameters.SORFactor
};
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
return apply.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
public JobHandle CalculateNormals(JobHandle inputDeps, float deltaTime)
{
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
var vorticity = new NormalsJob()
{
invMasses = solverImplementation.invMasses,
positions = solverImplementation.positions,
principalRadii = solverImplementation.principalRadii,
fluidMaterials = solverImplementation.fluidMaterials,
fluidMaterials2 = solverImplementation.fluidMaterials2,
fluidData = solverImplementation.fluidData,
fluidInterface = solverImplementation.fluidInterface,
velocities = solverImplementation.velocities,
angularVelocities = solverImplementation.angularVelocities,
vorticityAccelerations = solverImplementation.orientationDeltas.Reinterpret<float4>(),
vorticity = solverImplementation.restOrientations.Reinterpret<float4>(),
linearAccelerations = solverImplementation.positionDeltas,
linearFromAngular = solverImplementation.restPositions,
angularDiffusion = solverImplementation.anisotropies,
userData = solverImplementation.userData,
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
normals = solverImplementation.normals,
densityKernel = new Poly6Kernel(solverAbstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
gradKernel = new SpikyKernel(solverAbstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
solverParams = solverAbstraction.parameters,
batchData = batchData,
dt = deltaTime,
};
return vorticity.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
public JobHandle ViscosityAndVorticity(JobHandle inputDeps)
{
var eta = new ViscosityVorticityJob()
{
positions = solverImplementation.positions,
prevPositions = solverImplementation.prevPositions,
matchingRotations = solverImplementation.restPositions.Reinterpret<quaternion>(),
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
massCenters = solverImplementation.normals,
prevMassCenters = solverImplementation.renderablePositions,
fluidParams = solverImplementation.fluidMaterials,
deltas = solverImplementation.positionDeltas,
counts = solverImplementation.positionConstraintCounts,
batchData = batchData
};
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
return eta.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
public JobHandle AccumulateSmoothPositions(JobHandle inputDeps)
{
var accumulateSmooth = new AccumulateSmoothPositionsJob()
{
renderablePositions = solverImplementation.renderablePositions,
anisotropies = solverImplementation.anisotropies,
fluidMaterials = solverImplementation.fluidMaterials,
densityKernel = new Poly6Kernel(solverAbstraction.parameters.mode == Oni.SolverParameters.Mode.Mode2D),
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
batchData = batchData
};
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
return accumulateSmooth.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
public JobHandle AccumulateAnisotropy(JobHandle inputDeps)
{
var accumulateAnisotropy = new AccumulateAnisotropyJob()
{
renderablePositions = solverImplementation.renderablePositions,
anisotropies = solverImplementation.anisotropies,
pairs = ((BurstSolverImpl)constraints.solver).fluidInteractions,
batchData = batchData
};
int batchCount = batchData.isLast ? batchData.workItemCount : 1;
return accumulateAnisotropy.Schedule(batchData.workItemCount, batchCount, inputDeps);
}
[BurstCompile]
public struct UpdateDensitiesJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> positions;
[ReadOnly] public NativeArray<float4> prevPositions;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public NativeArray<FluidInteraction> pairs;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> moments;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> massCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> prevMassCenters;
[ReadOnly] public Poly6Kernel densityKernel;
[ReadOnly] public BatchData batchData;
[ReadOnly] public Oni.SolverParameters solverParams;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float restVolumeA = math.pow(principalRadii[pair.particleA].x * 2, 3 - (int)solverParams.mode);
float restVolumeB = math.pow(principalRadii[pair.particleB].x * 2, 3 - (int)solverParams.mode);
float gradA = restVolumeB * pair.avgGradient;
float gradB = restVolumeA * pair.avgGradient;
float vA = restVolumeB / restVolumeA;
float vB = restVolumeA / restVolumeB;
// accumulate pbf data (density, gradients):
fluidData[pair.particleA] += new float4(vA * pair.avgKernel, 0, gradA, gradA * gradA);
fluidData[pair.particleB] += new float4(vB * pair.avgKernel, 0, gradB, gradB * gradB);
// accumulate masses for COMs and moment matrices:
float wAvg = pair.avgKernel / ((densityKernel.W(0, fluidMaterials[pair.particleA].x) + densityKernel.W(0, fluidMaterials[pair.particleB].x)) * 0.5f);
massCenters[pair.particleA] += wAvg * new float4(positions[pair.particleB].xyz, 1) / positions[pair.particleB].w;
massCenters[pair.particleB] += wAvg * new float4(positions[pair.particleA].xyz, 1) / positions[pair.particleA].w;
prevMassCenters[pair.particleA] += wAvg * new float4(prevPositions[pair.particleB].xyz, 1) / positions[pair.particleB].w;
prevMassCenters[pair.particleB] += wAvg * new float4(prevPositions[pair.particleA].xyz, 1) / positions[pair.particleA].w;
moments[pair.particleA] += wAvg * (BurstMath.multrnsp4(positions[pair.particleB], prevPositions[pair.particleB]) + float4x4.identity * math.pow(principalRadii[pair.particleB].x, 2) * 0.001f) / positions[pair.particleB].w;
moments[pair.particleB] += wAvg * (BurstMath.multrnsp4(positions[pair.particleA], prevPositions[pair.particleA]) + float4x4.identity * math.pow(principalRadii[pair.particleA].x, 2) * 0.001f) / positions[pair.particleA].w;
}
}
}
[BurstCompile]
public struct ApplyDensityConstraintsJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public NativeArray<FluidInteraction> pairs;
[ReadOnly] public Poly6Kernel densityKernel;
[ReadOnly] public CohesionKernel cohesionKernel;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> positions;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[ReadOnly] public BatchData batchData;
[ReadOnly] public float sorFactor;
[ReadOnly] public Oni.SolverParameters solverParams;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float restVolumeA = math.pow(principalRadii[pair.particleA].x * 2, 3 - (int)solverParams.mode);
float restVolumeB = math.pow(principalRadii[pair.particleB].x * 2, 3 - (int)solverParams.mode);
float dist = math.length(positions[pair.particleA].xyz - positions[pair.particleB].xyz); // TODO: FIX! we cant read positions while we are writing to them.
// calculate tensile instability correction factor:
float cAvg = (cohesionKernel.W(dist, fluidMaterials[pair.particleA].x * 1.4f) + cohesionKernel.W(dist, fluidMaterials[pair.particleB].x * 1.4f)) * 0.5f;
float st = 0.2f * cAvg * (1 - math.saturate(math.abs(fluidMaterials[pair.particleA].y - fluidMaterials[pair.particleB].y))) * (fluidMaterials[pair.particleA].y + fluidMaterials[pair.particleB].y) * 0.5f;
float scorrA = -st / (positions[pair.particleA].w * fluidData[pair.particleA][3] + math.FLT_MIN_NORMAL);
float scorrB = -st / (positions[pair.particleB].w * fluidData[pair.particleB][3] + math.FLT_MIN_NORMAL);
// calculate position delta:
float4 delta = pair.gradient * pair.avgGradient * ((fluidData[pair.particleA][1] + scorrA) * restVolumeB + (fluidData[pair.particleB][1] + scorrB) * restVolumeA) * sorFactor;
delta.w = 0;
positions[pair.particleA] += delta * positions[pair.particleA].w;
positions[pair.particleB] -= delta * positions[pair.particleB].w;
}
}
}
[BurstCompile]
public struct ViscosityVorticityJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> positions;
[ReadOnly] public NativeArray<float4> prevPositions;
[ReadOnly] public NativeArray<quaternion> matchingRotations;
[ReadOnly] public NativeArray<float4> fluidParams;
[ReadOnly] public NativeArray<FluidInteraction> pairs;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> massCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> prevMassCenters;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> deltas;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<int> counts;
[ReadOnly] public BatchData batchData;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float visc = math.min(fluidParams[pair.particleA].z, fluidParams[pair.particleB].z);
// viscosity:
float4 goalA = new float4(massCenters[pair.particleB].xyz + math.rotate(matchingRotations[pair.particleB], (prevPositions[pair.particleA] - prevMassCenters[pair.particleB]).xyz), 0);
float4 goalB = new float4(massCenters[pair.particleA].xyz + math.rotate(matchingRotations[pair.particleA], (prevPositions[pair.particleB] - prevMassCenters[pair.particleA]).xyz), 0);
deltas[pair.particleA] += (goalA - positions[pair.particleA]) * visc;
deltas[pair.particleB] += (goalB - positions[pair.particleB]) * visc;
counts[pair.particleA]++;
counts[pair.particleB]++;
}
}
}
[BurstCompile]
public struct NormalsJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float> invMasses;
[ReadOnly] public NativeArray<float4> velocities;
[ReadOnly] public NativeArray<float4> angularVelocities;
[ReadOnly] public NativeArray<float4> positions;
[ReadOnly] public NativeArray<float4> vorticity;
[ReadOnly] public NativeArray<float4> principalRadii;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public NativeArray<float4> fluidMaterials2;
[ReadOnly] public NativeArray<float4> fluidInterface;
[ReadOnly] public NativeArray<FluidInteraction> pairs;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> fluidData;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> userData;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4> normals;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> linearAccelerations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> vorticityAccelerations;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4> linearFromAngular;
[NativeDisableContainerSafetyRestriction] [NativeDisableParallelForRestriction] public NativeArray<float4x4> angularDiffusion;
[ReadOnly] public Poly6Kernel densityKernel;
[ReadOnly] public SpikyKernel gradKernel;
[ReadOnly] public BatchData batchData;
[ReadOnly] public Oni.SolverParameters solverParams;
[ReadOnly] public float dt;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float restVolumeA = math.pow(principalRadii[pair.particleA].x * 2, 3 - (int)solverParams.mode);
float restVolumeB = math.pow(principalRadii[pair.particleB].x * 2, 3 - (int)solverParams.mode);
float invDensityA = invMasses[pair.particleA] / fluidData[pair.particleA].x;
float invDensityB = invMasses[pair.particleB] / fluidData[pair.particleB].x;
float3 relVel = velocities[pair.particleA].xyz - velocities[pair.particleB].xyz;
float3 relAng = angularVelocities[pair.particleA].xyz - angularVelocities[pair.particleB].xyz;
float3 relVort = vorticity[pair.particleA].xyz - vorticity[pair.particleB].xyz;
float4 d = new float4((positions[pair.particleA] - positions[pair.particleB]).xyz,0);
float dist = math.length(d);
float avgGrad = (gradKernel.W(dist, fluidMaterials[pair.particleA].x) +
gradKernel.W(dist, fluidMaterials[pair.particleB].x)) * 0.5f;
float avgKern = (densityKernel.W(dist, fluidMaterials[pair.particleA].x) +
densityKernel.W(dist, fluidMaterials[pair.particleB].x)) * 0.5f;
float avgNorm = (densityKernel.W(0, fluidMaterials[pair.particleA].x) +
densityKernel.W(0, fluidMaterials[pair.particleB].x)) * 0.5f;
// property diffusion:
float diffusionSpeed = (fluidInterface[pair.particleA].w + fluidInterface[pair.particleB].w) * avgKern * dt;
float4 userDelta = (userData[pair.particleB] - userData[pair.particleA]) * solverParams.diffusionMask * diffusionSpeed;
userData[pair.particleA] += restVolumeB / restVolumeA * userDelta;
userData[pair.particleB] -= restVolumeA / restVolumeB * userDelta;
// calculate color field normal:
float4 normGrad = d / (dist + BurstMath.epsilon);
float4 vgrad = normGrad * avgGrad;
float radius = (fluidMaterials[pair.particleA].x + fluidMaterials[pair.particleB].x) * 0.5f;
normals[pair.particleA] += vgrad * radius * restVolumeB;
normals[pair.particleB] -= vgrad * radius * restVolumeA;
// measure relative velocity for foam generation:
float4 dataA = fluidData[pair.particleA];
float4 dataB = fluidData[pair.particleB];
float relVelMag = math.length(relVel) + BurstMath.epsilon;
float avgVelDiffKernel = 1 - math.min(1, dist / (radius + BurstMath.epsilon));
float rv = relVelMag * (1 - math.dot(relVel / relVelMag, normGrad.xyz)) * avgVelDiffKernel;
dataA.z += rv;
dataB.z += rv;
fluidData[pair.particleA] = dataA;
fluidData[pair.particleB] = dataB;
// micropolar: curl of linear/angular velocity:
float3 velCross = math.cross(relVel, vgrad.xyz);
float3 vortCross = math.cross(relVort, vgrad.xyz);
linearAccelerations[pair.particleA] += new float4(vortCross / invMasses[pair.particleB] * invDensityA, 0);
linearAccelerations[pair.particleB] += new float4(vortCross / invMasses[pair.particleA] * invDensityB, 0);
vorticityAccelerations[pair.particleA] += new float4(velCross / invMasses[pair.particleB] * invDensityA, 0);
vorticityAccelerations[pair.particleB] += new float4(velCross / invMasses[pair.particleA] * invDensityB, 0);
// angular diffusion:
float4x4 diffA = angularDiffusion[pair.particleA];
float4x4 diffB = angularDiffusion[pair.particleB];
diffA.c0 += new float4(relAng * avgKern / invMasses[pair.particleB] * invDensityA, 0);
diffB.c0 -= new float4(relAng * avgKern / invMasses[pair.particleA] * invDensityB, 0);
diffA.c1 += new float4(relVort * avgKern / invMasses[pair.particleB] * invDensityA, 0);
diffB.c1 -= new float4(relVort * avgKern / invMasses[pair.particleA] * invDensityB, 0);
angularDiffusion[pair.particleA] = diffA;
angularDiffusion[pair.particleB] = diffB;
// linear velocity due to baroclinity:
linearFromAngular[pair.particleA] += new float4(math.cross(angularVelocities[pair.particleB].xyz, d.xyz) * avgKern / avgNorm, 0);
linearFromAngular[pair.particleB] -= new float4(math.cross(angularVelocities[pair.particleA].xyz, d.xyz) * avgKern / avgNorm, 0);
}
}
}
[BurstCompile]
public struct AccumulateSmoothPositionsJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> renderablePositions;
[ReadOnly] public NativeArray<float4> fluidMaterials;
[ReadOnly] public Poly6Kernel densityKernel;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4x4> anisotropies;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<FluidInteraction> pairs;
[ReadOnly] public BatchData batchData;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float distance = math.length((renderablePositions[pair.particleA] - renderablePositions[pair.particleB]).xyz);
pair.avgKernel = (densityKernel.W(distance, fluidMaterials[pair.particleA].x) +
densityKernel.W(distance, fluidMaterials[pair.particleB].x)) * 0.5f;
var A = anisotropies[pair.particleA];
var B = anisotropies[pair.particleB];
A.c3 += new float4(renderablePositions[pair.particleB].xyz,1) * pair.avgKernel;
B.c3 += new float4(renderablePositions[pair.particleA].xyz,1) * pair.avgKernel;
anisotropies[pair.particleA] = A;
anisotropies[pair.particleB] = B;
pairs[i] = pair;
}
}
}
[BurstCompile]
public struct AccumulateAnisotropyJob : IJobParallelFor
{
[ReadOnly] public NativeArray<float4> renderablePositions;
[ReadOnly] public NativeArray<FluidInteraction> pairs;
[NativeDisableContainerSafetyRestriction][NativeDisableParallelForRestriction] public NativeArray<float4x4> anisotropies;
[ReadOnly] public BatchData batchData;
public void Execute(int workItemIndex)
{
int start, end;
batchData.GetConstraintRange(workItemIndex, out start, out end);
for (int i = start; i < end; ++i)
{
var pair = pairs[i];
float4 distanceA = (renderablePositions[pair.particleB] - anisotropies[pair.particleA].c3) * pair.avgKernel;
float4 distanceB = (renderablePositions[pair.particleA] - anisotropies[pair.particleB].c3) * pair.avgKernel;
anisotropies[pair.particleA] += BurstMath.multrnsp4(distanceA,distanceA);
anisotropies[pair.particleB] += BurstMath.multrnsp4(distanceB,distanceB);
}
}
}
}
}
#endif

View File

@@ -1,11 +0,0 @@
fileFormatVersion: 2
guid: 8292e6ef6129f47abaee4fb2cb49055e
MonoImporter:
externalObjects: {}
serializedVersion: 2
defaultReferences: []
executionOrder: 0
icon: {instanceID: 0}
userData:
assetBundleName:
assetBundleVariant:

View File

@@ -1,8 +0,0 @@
fileFormatVersion: 2
guid: 5a128a7c745c84794a944362f49011fc
folderAsset: yes
DefaultImporter:
externalObjects: {}
userData:
assetBundleName:
assetBundleVariant:

View File

@@ -1,15 +0,0 @@
#if (OBI_BURST && OBI_MATHEMATICS && OBI_COLLECTIONS)
using System;
using Unity.Mathematics;
namespace Obi
{
public struct CohesionKernel
{
public float W(float r, float h)
{
return math.cos(math.min(r, h) * 3 * math.PI / (2 * h));
}
}
}
#endif

View File

@@ -1,11 +0,0 @@
fileFormatVersion: 2
guid: d2bcf975ae59e4fcf9412ec328a3e2b9
MonoImporter:
externalObjects: {}
serializedVersion: 2
defaultReferences: []
executionOrder: 0
icon: {instanceID: 0}
userData:
assetBundleName:
assetBundleVariant:

View File

@@ -1,36 +0,0 @@
#if (OBI_BURST && OBI_MATHEMATICS && OBI_COLLECTIONS)
using System;
using Unity.Mathematics;
namespace Obi
{
public struct Poly6Kernel
{
public float norm;
public bool norm2D;
public Poly6Kernel(bool norm2D)
{
this.norm2D = norm2D;
if (norm2D)
norm = 4.0f / math.PI;
else
norm = 315.0f / (64.0f * math.PI);
}
public float W(float r, float h)
{
float h2 = h * h;
float h4 = h2 * h2;
float h8 = h4 * h4;
float rl = math.min(r, h);
float hr = h2 - rl * rl;
if (norm2D)
return norm / h8 * hr * hr * hr;
return norm / (h8 * h) * hr * hr * hr;
}
}
}
#endif

View File

@@ -1,11 +0,0 @@
fileFormatVersion: 2
guid: 7a3990134524143ac852b488554f6d4e
MonoImporter:
externalObjects: {}
serializedVersion: 2
defaultReferences: []
executionOrder: 0
icon: {instanceID: 0}
userData:
assetBundleName:
assetBundleVariant:

View File

@@ -1,35 +0,0 @@
#if (OBI_BURST && OBI_MATHEMATICS && OBI_COLLECTIONS)
using System;
using Unity.Mathematics;
namespace Obi
{
public struct SpikyKernel
{
public float norm;
public bool norm2D;
public SpikyKernel(bool norm2D)
{
this.norm2D = norm2D;
if (norm2D)
norm = -30.0f / math.PI;
else
norm = -45.0f / math.PI;
}
public float W(float r, float h)
{
float h2 = h * h;
float h4 = h2 * h2;
float rl = math.min(r, h);
float hr = h - rl;
if (norm2D)
return norm / (h4 * h) * hr * hr;
return norm / (h4 * h2) * hr * hr;
}
}
}
#endif

View File

@@ -1,11 +0,0 @@
fileFormatVersion: 2
guid: a407989bfa0664e9ab75773d1808f549
MonoImporter:
externalObjects: {}
serializedVersion: 2
defaultReferences: []
executionOrder: 0
icon: {instanceID: 0}
userData:
assetBundleName:
assetBundleVariant: