| 44 | 
  | 
 * | 
| 45 | 
  | 
 *  Created by Charles F. Vardeman II on 29 July 2008. | 
| 46 | 
  | 
 *  @author  Charles F. Vardeman II | 
| 47 | 
< | 
 *  @version $Id: Triangle.cpp,v 1.1 2008-11-14 15:44:34 chuckv Exp $ | 
| 47 | 
> | 
 *  @version $Id: Triangle.cpp,v 1.3 2009-05-13 22:27:29 gezelter Exp $ | 
| 48 | 
  | 
 * | 
| 49 | 
  | 
 */ | 
| 50 | 
  | 
 | 
| 53 | 
  | 
using namespace oopse; | 
| 54 | 
  | 
 | 
| 55 | 
  | 
 | 
| 56 | 
< | 
Triangle::Triangle() : HaveNormal_(false), HaveCentroid_(false),HaveArea_(false), area_(0.0), normal_(V3Zero),  | 
| 57 | 
< | 
                       centroid_(V3Zero),facetVelocity_(V3Zero), mass_(0.0), | 
| 58 | 
< | 
                       a_(V3Zero),b_(V3Zero),c_(V3Zero){ | 
| 56 | 
> | 
Triangle::Triangle() : HaveNormal_(false), HaveCentroid_(false), | 
| 57 | 
> | 
                       HaveArea_(false), area_(0.0), normal_(V3Zero),  | 
| 58 | 
> | 
                       centroid_(V3Zero), facetVelocity_(V3Zero), mass_(0.0), | 
| 59 | 
> | 
                       a_(V3Zero), b_(V3Zero), c_(V3Zero){ | 
| 60 | 
  | 
} | 
| 61 | 
  | 
 | 
| 62 | 
  | 
void Triangle::addVertices(Vector3d P1, Vector3d P2, Vector3d P3){ | 
| 67 | 
  | 
  // Compute some quantites like a,b,c | 
| 68 | 
  | 
  a_ = P1-P2; | 
| 69 | 
  | 
  b_ = P1-P3; | 
| 70 | 
< | 
  c_ = P1-P3;   | 
| 70 | 
> | 
  c_ = P2-P3;   | 
| 71 | 
  | 
} | 
| 72 | 
  | 
 | 
| 73 | 
  | 
 | 
| 88 | 
  | 
  centroid_ = (vertices_[0] + vertices_[1] + vertices_[2])/3.0; | 
| 89 | 
  | 
  return centroid_; | 
| 90 | 
  | 
} | 
| 91 | 
+ | 
 | 
| 92 | 
+ | 
 | 
| 93 | 
+ | 
Mat3x3d Triangle::computeHydrodynamicTensor(RealType viscosity) { | 
| 94 | 
+ | 
   | 
| 95 | 
+ | 
  Vector3d u0 = -a_; | 
| 96 | 
+ | 
  Vector3d v0 = centroid_ - vertices_[0]; | 
| 97 | 
+ | 
  RealType s0 = 0.5*cross(u0,v0).length(); | 
| 98 | 
+ | 
   | 
| 99 | 
+ | 
  Vector3d u1 = -c_; | 
| 100 | 
+ | 
  Vector3d v1 = centroid_ - vertices_[1]; | 
| 101 | 
+ | 
  RealType s1 = 0.5*cross(u1,v1).length(); | 
| 102 | 
+ | 
   | 
| 103 | 
+ | 
  Vector3d u2 = b_; | 
| 104 | 
+ | 
  Vector3d v2 = centroid_ - vertices_[2]; | 
| 105 | 
+ | 
  RealType s2 = 0.5*cross(u2,v2).length(); | 
| 106 | 
+ | 
   | 
| 107 | 
+ | 
  Mat3x3d H; | 
| 108 | 
+ | 
  H = hydro_tensor(centroid_,centroid_,vertices_[1],vertices_[0],s0,viscosity)+ | 
| 109 | 
+ | 
    hydro_tensor(centroid_,centroid_,vertices_[1],vertices_[2],s1,viscosity)+ | 
| 110 | 
+ | 
    hydro_tensor(centroid_,centroid_,vertices_[2],vertices_[0],s2,viscosity); | 
| 111 | 
+ | 
 | 
| 112 | 
+ | 
  return H.inverse(); | 
| 113 | 
+ | 
} | 
| 114 | 
+ | 
 | 
| 115 | 
+ | 
Mat3x3d Triangle::hydro_tensor( | 
| 116 | 
+ | 
                               const Vector3d& ri, | 
| 117 | 
+ | 
                               const Vector3d& rj0, | 
| 118 | 
+ | 
                               const Vector3d& rj1, | 
| 119 | 
+ | 
                               const Vector3d& rj2, | 
| 120 | 
+ | 
                               RealType s, RealType viscosity){ | 
| 121 | 
+ | 
   | 
| 122 | 
+ | 
  Vector3d v2 = (rj0 + rj1 + rj2)/3.0;  // sub-centroid | 
| 123 | 
+ | 
  Vector3d dr = ri - v2;                // real centroid to sub-centroid | 
| 124 | 
+ | 
  RealType l2 = 1.0/dr.lengthSquare();   | 
| 125 | 
+ | 
  | 
| 126 | 
+ | 
  Mat3x3d G; | 
| 127 | 
+ | 
  G = (SquareMatrix3<RealType>::identity() + outProduct(dr,dr)*l2)*sqrt(l2); | 
| 128 | 
+ | 
 | 
| 129 | 
+ | 
  G *= 0.125/3.14159285358979; | 
| 130 | 
+ | 
  G *= s/viscosity; | 
| 131 | 
+ | 
  return G; | 
| 132 | 
+ | 
} |