00001 #ifndef OPENTISSUE_DYNAMICS_SPH_SOLVERS_SPH_VISCOSITY_H 00002 #define OPENTISSUE_DYNAMICS_SPH_SOLVERS_SPH_VISCOSITY_H 00003 // 00004 // OpenTissue Template Library 00005 // - A generic toolbox for physics-based modeling and simulation. 00006 // Copyright (C) 2008 Department of Computer Science, University of Copenhagen. 00007 // 00008 // OTTL is licensed under zlib: http://opensource.org/licenses/zlib-license.php 00009 // 00010 #include <OpenTissue/configuration.h> 00011 00012 #include <OpenTissue/dynamics/sph/sph_solver.h> 00013 00014 namespace OpenTissue 00015 { 00016 namespace sph 00017 { 00018 00022 template< typename Types, typename KernelPolicy > 00023 class ViscosityForce : Solver< Types, typename Types::vector > 00024 { 00025 public: 00026 typedef Solver<Types, typename Types::vector> base_type; 00027 typedef KernelPolicy smoothing_kernel; 00028 typedef typename base_type::value value; 00029 typedef typename Types::real_type real_type; 00030 typedef typename Types::particle particle; 00031 typedef typename Types::particle_cptr_container::const_iterator particle_cptr_container_citerator; 00032 00033 public: 00037 ViscosityForce(const real_type& viscosity = 0) : base_type() 00038 { 00039 setViscosity(viscosity); 00040 } 00041 00045 ~ViscosityForce() 00046 { 00047 } 00048 00052 ViscosityForce& operator=(const ViscosityForce& rhs) 00053 { 00054 m_m = rhs.m_m; 00055 return *this; 00056 } 00057 00058 public: 00059 void setViscosity(const real_type& coefficient) 00060 { 00061 m_m = coefficient; 00062 } 00063 00069 virtual value apply(const particle& par, particle_cptr_container_citerator begin, particle_cptr_container_citerator end) const 00070 { 00071 #if 1 00072 value res(0); // reset vector; 00073 for (particle_cptr_container_citerator p_ = begin; p_ != end; ++p_) { 00074 const particle* p = *p_; 00075 if (&par == p) continue; // all sph forces apply i != j 00076 res += p->mass()*((p->velocity()-par.velocity())/p->density())*m_W.laplacian(par.position()-p->position()); 00077 } 00078 return value(m_m*res); 00079 #else 00080 value res(0); // reset vector; 00081 for (particle_cptr_container_citerator p_ = begin; p_ != end; ++p_) { 00082 const particle* p = *p_; 00083 if (&par == p) continue; // all sph forces apply i != j 00084 res += p->mass()*(p->velocity()-par.velocity())*m_W.laplacian(par.position()-p->position()); 00085 } 00086 return value(m_m*res/par.density()); 00087 #endif 00088 } 00089 00090 00096 value apply(const particle& par, const particle& p) const 00097 { 00098 return value(m_m*(p.mass()*((p.velocity()-par.velocity())/p.density())*m_W.laplacian(par.position()-p.position()))); 00099 } 00100 00101 private: 00102 const smoothing_kernel m_W; 00103 real_type m_m; 00104 00105 }; // End class ViscosityForce 00106 00107 } // namespace sph 00108 } // namespace OpenTissue 00109 00110 // OPENTISSUE_DYNAMICS_SPH_SOLVERS_SPH_VISCOSITY_H 00111 #endif
1.7.1