[osg-users] Fluid Friction

[Max Bandazian]

I've been attempting to use osg particle systems, and I believe I've
encountered either a bug or some undocumented behavior with fluid friction.

Using, for instance, a particle system configured to look smoke-plume-like,
and using a ModularProgram with a FluidFrictionOperator, I am finding that
sometimes the plume appears normally, and sometimes the smoke particles fly
off into space at high velocity in a spiral pattern. I tried a variety of
things to track down the problem; here's what I did, with notes:

- Removed the friction operator:
  The smoke displayed correctly, but of course frictionless.

- Replaced the ModularProgram with a FluidProgram:
  The smoke displayed correctly, with friction.

- Replaced the FluidFrictionOperator::operate code with the calculation from
FluidProgram:
  The smoke displayed correctly, with friction. The FluidProgram's
calculation appears totally different from
  the FluidFrictionOperator's calculation, and uses different variable names
(coeffA/coeffB vs. viscosityCoefficient/densityCoefficient)
  for what looks like the same thing. It's hard to see why.

- Replaced the FluidFrictionOperator::operate code with a quadratic drag
calculation from a physics book:
  The smoke displayed correctly, with friction.

So I am wondering if anyone can enlighten me as to what is going on here.
Why are the FluidProgram and FluidFrictionOperator calculations so
different, what's with the cryptic variable names in FluidFrictionOperator,
and why would one calculation work but not the other?

I'm including the friction calculation I settled on, in case anyone is
curious or would like to point out the many mistakes I'm sure I made.

Thanks for any help.
Max Bandazian

   osg::Vec3 computeQuadraticDrag(osgParticle::Particle* particle, double
dt,
      float viscosityCoefficient, float densityCoefficient, osg::Vec3 wind)
   {
      static const float four_over_three = 4.0/3.0;

      float radius = particle->getRadius();
      float area = osg::PI*radius*radius;
      float volume = area*radius*four_over_three;

      osg::Vec3 vRelativeToFluid(particle->getVelocity() - wind);
      float speedRelativeToFluid = vRelativeToFluid.length();

      osg::Vec3 vDirection(particle->getVelocity());
      vDirection.normalize();

      float dragCoefficient = 6 * osg::PI * viscosityCoefficient * radius;
// Stokes

      osg::Vec3 dragForce = vDirection * -0.5 * densityCoefficient *
         speedRelativeToFluid * speedRelativeToFluid * area;

      return dragForce;
   }

   osg::Vec3 computeWindForce(osgParticle::Particle* particle, double dt,
      float viscosityCoefficient, float densityCoefficient, osg::Vec3 wind)
   {
      static const float four_over_three = 4.0/3.0;

      float radius = particle->getRadius();
      float area = osg::PI*radius*radius;
      float volume = area*radius*four_over_three;

      // Find wind vector minus projection of particle velocity
      osg::Vec3 windNorm(wind);
      windNorm.normalize();
      float windDot = windNorm * particle->getVelocity();
      osg::Vec3 relativeWind = wind * (1 - windDot);

      osg::Vec3 windForce = (relativeWind - particle->getVelocity())
         * densityCoefficient * volume / dt;

      return windForce;
   }

   void MyFluidFrictionOperator::operate(osgParticle::Particle* particle,
double dt)
   {
      osg::Vec3 dragForce = computeQuadraticDrag(particle, dt,
viscosityCoefficient,
         densityCoefficient, _wind);
      osg::Vec3 windForce = computeWindForce(particle, dt,
viscosityCoefficient,
         densityCoefficient, _wind);

      osg::Vec3 acceleration = (dragForce + windForce) *
particle->getMassInv();
      osg::Vec3 velocity = acceleration * dt;

      particle->addVelocity(velocity);
   }

_______________________________________________
osg-users mailing list
osg-users@lists.openscenegraph.org
http://lists.openscenegraph.org/listinfo.cgi/osg-users-openscenegraph.org