ESyS-Particle
4.0.1
|
00001 00002 // // 00003 // Copyright (c) 2003-2011 by The University of Queensland // 00004 // Earth Systems Science Computational Centre (ESSCC) // 00005 // http://www.uq.edu.au/esscc // 00006 // // 00007 // Primary Business: Brisbane, Queensland, Australia // 00008 // Licensed under the Open Software License version 3.0 // 00009 // http://www.opensource.org/licenses/osl-3.0.php // 00010 // // 00012 00013 #include "Foundation/console.h" 00014 00022 template<class ParticleType,class IType> 00023 TriMesh_PIS_NE<ParticleType,IType>::TriMesh_PIS_NE(TriMesh* mesh_p,ParallelParticleArray<ParticleType>* ppa_p,typename IType::ParameterType param) 00024 :TriMesh_PIS<ParticleType>(mesh_p,ppa_p) 00025 { 00026 m_param=param; 00027 this->m_update_timestamp=0; 00028 } 00029 00033 template<class ParticleType,class IType> 00034 TriMesh_PIS_NE<ParticleType,IType>::~TriMesh_PIS_NE() 00035 {} 00036 00046 template <class ParticleType,class IType> 00047 bool TriMesh_PIS_NE<ParticleType,IType>::isIn(const std::vector<int>& v) 00048 { 00049 bool res=false; 00050 00051 if(v.size()<3){ 00052 res=false; 00053 } else { 00054 switch (v[2]){ 00055 case 0: res=m_tri_int_set.find(make_pair(v[0],v[1]))!=m_tri_int_set.end(); break; 00056 case 1: res=m_edge_int_set.find(make_pair(v[0],v[1]))!=m_edge_int_set.end(); break; 00057 case 2: res=m_corner_int_set.find(make_pair(v[0],v[1]))!=m_corner_int_set.end(); break; 00058 default: console.Error() << "wrong value in argument of TriMesh_PIS::isIn !!\n"; break; 00059 } 00060 } 00061 00062 return res; 00063 } 00064 00068 template<class ParticleType,class IType> 00069 void TriMesh_PIS_NE<ParticleType,IType>::calcForces() 00070 { 00071 console.XDebug() << "TriMesh_PIS_NE calculating " << m_triangle_interactions.size() << " triangle forces , " 00072 << m_edge_interactions.size() << "edge forces and " 00073 << m_corner_interactions.size() << "corner forces\n"; 00074 00075 // calculate forces for triangle interactions 00076 for(typename std::vector<typename IType::TriIntType>::iterator tri_iter=m_triangle_interactions.begin(); 00077 tri_iter!=m_triangle_interactions.end(); 00078 tri_iter++){ 00079 tri_iter->calcForces(); 00080 } 00081 // calculate forces for edge interactions 00082 for(typename std::vector<typename IType::EdgeIntType>::iterator edge_iter=m_edge_interactions.begin(); 00083 edge_iter!=m_edge_interactions.end(); 00084 edge_iter++){ 00085 edge_iter->calcForces(); 00086 } 00087 // calculate forces for corner interactions 00088 for(typename std::vector<typename IType::CornerIntType>::iterator corner_iter=m_corner_interactions.begin(); 00089 corner_iter!=m_corner_interactions.end(); 00090 corner_iter++){ 00091 corner_iter->calcForces(); 00092 } 00093 } 00094 00098 template<class ParticleType,class IType> 00099 bool TriMesh_PIS_NE<ParticleType,IType>::update() 00100 { 00101 console.XDebug() << "TriMesh_PIS_NE::update\n"; 00102 bool res=false; 00103 //int count_edge=0; 00104 //int count_tri=0; 00105 00106 if(this->m_update_timestamp != this->m_ppa->getTimeStamp()){// m_ppa rebuild since last update 00107 console.XDebug() << "TriMesh_PIS_NE doing update\n"; 00108 // clean out old interactions 00109 m_triangle_interactions.clear(); 00110 m_edge_interactions.clear(); 00111 m_corner_interactions.clear(); 00112 m_tri_int_set.clear(); 00113 m_edge_int_set.clear(); 00114 m_corner_int_set.clear(); 00115 // -- get triangle interactions -- 00116 // for all triangles 00117 for( 00118 TriMesh::triangle_iterator tri_iter = this->m_mesh->triangles_begin(); 00119 tri_iter != this->m_mesh->triangles_end(); 00120 tri_iter++){ 00121 // get list of particles near this triangle 00122 typename ParallelParticleArray<ParticleType>::ParticleListHandle plh= 00123 ((ParallelParticleArray<ParticleType>*)this->m_ppa)->getParticlesNearTriangle(*tri_iter); 00124 // for all particles found 00125 for( 00126 typename ParallelParticleArray<ParticleType>::ParticleListIterator p_iter=plh->begin(); 00127 p_iter!=plh->end(); 00128 p_iter++ 00129 ){ 00130 bool iflag = this->m_ppa->isInInner((*p_iter)->getPos()); 00131 m_triangle_interactions.push_back(typename IType::TriIntType(*p_iter,&(*tri_iter),m_param,iflag)); 00132 m_tri_int_set.insert(make_pair(tri_iter->getID(),(*p_iter)->getID())); 00133 } 00134 } 00135 // -- get edge interactions 00136 // for all edges 00137 for( 00138 TriMesh::edge_iterator ed_iter = this->m_mesh->edges_begin(); 00139 ed_iter != this->m_mesh->edges_end(); 00140 ed_iter++ 00141 ){ 00142 typename ParallelParticleArray<ParticleType>::ParticleListHandle plh = 00143 ((ParallelParticleArray<ParticleType>*)(this->m_ppa))->getParticlesNearEdge(&(*ed_iter)); 00144 for ( 00145 typename ParallelParticleArray<ParticleType>::ParticleListIterator p_iter=plh->begin(); 00146 p_iter!=plh->end(); 00147 p_iter++ 00148 ){ 00149 bool iflag = this->m_ppa->isInInner((*p_iter)->getPos()); 00150 m_edge_interactions.push_back(typename IType::EdgeIntType(*p_iter,&(*ed_iter),m_param,iflag)); 00151 //m_particle_id_set.insert((*p_iter)->getID()); 00152 } 00153 } 00154 // --- get corner interactions 00155 for ( 00156 TriMesh::corner_iterator co_iter = this->m_mesh->corners_begin(); 00157 co_iter != this->m_mesh->corners_end(); 00158 co_iter++ 00159 ){ 00160 typename ParallelParticleArray<ParticleType>::ParticleListHandle plh= 00161 ((ParallelParticleArray<ParticleType>*)(this->m_ppa))->getParticlesNearPoint(co_iter->getPos()); 00162 for ( 00163 typename ParallelParticleArray<ParticleType>::ParticleListIterator p_iter=plh->begin(); 00164 p_iter!=plh->end(); 00165 p_iter++ 00166 ){ 00167 bool iflag = this->m_ppa->isInInner((*p_iter)->getPos()); 00168 m_corner_interactions.push_back(typename IType::CornerIntType(*p_iter,&(*co_iter),m_param,iflag)); 00169 //m_particle_id_set.insert((*p_iter)->getID()); 00170 } 00171 } 00172 // set timestamp 00173 this->m_update_timestamp = this->m_ppa->getTimeStamp(); 00174 } 00175 console.XDebug() << "end ElasticTriMeshIG<T>::Update\n"; 00176 00177 return res; 00178 } 00179 00180