ESyS-Particle  4.0.1
SoftBWallInteraction.hpp
1 
2 // //
3 // Copyright (c) 2003-2011 by The University of Queensland //
4 // Earth Systems Science Computational Centre (ESSCC) //
5 // http://www.uq.edu.au/esscc //
6 // //
7 // Primary Business: Brisbane, Queensland, Australia //
8 // Licensed under the Open Software License version 3.0 //
9 // http://www.opensource.org/licenses/osl-3.0.php //
10 // //
12 
13 
14 #ifndef MODEL_SOFTBWALLINTERACTION_HPP
15 #define MODEL_SOFTBWALLINTERACTION_HPP
16 
17 #include "SoftBWallInteraction.h"
18 
19 template <class T>
20 CSoftBondedWallInteraction<T>::CSoftBondedWallInteraction(T* p,CWall* w,double normalK,double shearK,bool scaling,bool iflag):
21  AWallInteraction<T>(p,w,iflag)
22 {
23  double scale;
24 
25  if (scaling) {
26  // scale stiffness to particle cross section
27  if(CParticle::getDo2dCalculations()){ // 2D
28 // scale=2.0*this->m_p->getRad();
29  scale=1.0;
30  } else { // 3D
31 // scale=3.1415926536*this->m_p->getRad()*this->m_p->getRad();
32  scale=3.1415926536*this->m_p->getRad();
33  }
34  }
35  else {
36  scale = 1.0;
37  }
38 
39  m_normalK=normalK*scale;
40  m_shearK=shearK*scale;
41 }
42 
46 template <class T>
48 {
49  const Vec3 &n = this->m_wall->getNormal();
50  const Vec3 relDisplacement =
51  (
52  this->m_p->getTotalDisplacement()
53  -
54  this->m_wall->getTotalDisplacement()
55  );
56 
57  const double normalDist = dot(relDisplacement, n)/(n.norm());
58  const Vec3 normalForce = ((m_normalK*normalDist)/n.norm())*n;
59 
60  /*
61  * Tangent displacement vector is relDisplacement point
62  * projected onto plane.
63  */
64  const Vec3 tangentDisplacement =
65  (relDisplacement - ((normalDist/(n.norm()))*n));
66  const Vec3 tangentForce = m_shearK*tangentDisplacement;
67 
68  const Vec3 totalForce = normalForce + tangentForce;
69  this->m_p->applyForce(-1.0*totalForce,this->m_p->getPos());
70  if(this->m_inner_flag) this->m_wall->addForce(totalForce);
71 }
72 
73 #endif