1 #include <osgParticle/FluidProgram> 2 FluidProgram()3osgParticle::FluidProgram::FluidProgram(): 4 Program() 5 { 6 setFluidToAir(); 7 } 8 FluidProgram(const FluidProgram & copy,const osg::CopyOp & copyop)9osgParticle::FluidProgram::FluidProgram(const FluidProgram& copy, const osg::CopyOp& copyop): 10 Program(copy, copyop), 11 _acceleration(copy._acceleration), 12 _viscosity(copy._viscosity), 13 _density(copy._density), 14 _wind(copy._wind), 15 _viscosityCoefficient(copy._viscosityCoefficient), 16 _densityCoefficient(copy._densityCoefficient) 17 { 18 } 19 execute(double dt)20void osgParticle::FluidProgram::execute(double dt) 21 { 22 const float four_over_three = 4.0f/3.0f; 23 ParticleSystem* ps = getParticleSystem(); 24 int n = ps->numParticles(); 25 for (int i=0; i<n; ++i) 26 { 27 Particle* particle = ps->getParticle(i); 28 if (particle->isAlive()) 29 { 30 float radius = particle->getRadius(); 31 float Area = osg::PI*radius*radius; 32 float Volume = Area*radius*four_over_three; 33 34 // compute force due to gravity + boyancy of displacing the fluid that the particle is emersed in. 35 osg::Vec3 accel_gravity = _acceleration * ((particle->getMass() - _density*Volume) * particle->getMassInv()); 36 37 // compute force due to friction 38 osg::Vec3 relative_wind = particle->getVelocity()-_wind; 39 osg::Vec3 wind_force = - relative_wind * Area * (_viscosityCoefficient + _densityCoefficient*relative_wind.length()); 40 osg::Vec3 wind_accel = wind_force * particle->getMassInv(); 41 42 double compenstated_dt = dt; 43 if (relative_wind.length2() < dt*dt*wind_accel.length2()) 44 { 45 // OSG_NOTICE<<"** Could be critical: dt="<<dt<<" critical_dt="<<sqrtf(critical_dt2)<<std::endl; 46 double critical_dt2 = relative_wind.length2()/wind_accel.length2(); 47 compenstated_dt = sqrtf(critical_dt2)*0.8f; 48 } 49 50 particle->addVelocity(accel_gravity*dt + wind_accel*compenstated_dt); 51 52 53 } 54 } 55 } 56