diff --git a/contrib/arc/Classes/NodeEnrichedFS.cpp b/contrib/arc/Classes/NodeEnrichedFS.cpp new file mode 100644 index 0000000000000000000000000000000000000000..44bd73527f0f57753d53951b290459c4baec5c7e --- /dev/null +++ b/contrib/arc/Classes/NodeEnrichedFS.cpp @@ -0,0 +1,254 @@ +// +// Description : From function space to node enriched composant function space +// +// +// Author: <Eric Bechet>::<Boris Sedji>, 01/2010 +// +// Copyright: See COPYING file that comes with this distribution +// +// + + +#include "NodeEnrichedFS.h" + +template <class T> void NodeEnrichedFS<T>::f(MVertex *ver, std::vector<ValType> &vals) +{ + std::vector<ValType> valsd; + // get the support value at the vertex + _SupportFS->f(ver,valsd); + int normaldofs=valsd.size(); + int curpos = vals.size(); + + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(ver->getNum()); + if (it!=_TagEnrichedVertex->end()) // if node enriched + { + double func = (*_funcEnrichment)(ver->x(),ver->y(),ver->z()); + int nbdofs = normaldofs + _EnrichComp->size(); + vals.reserve(curpos+nbdofs); + for (int i=0;i<normaldofs;i++) vals.push_back(valsd[i]); + for (int i=0;i<_EnrichComp->size();i++) vals.push_back(valsd[_EnrichComp->at(i)]*func); + } + else + { + int nbdofs = normaldofs; + vals.reserve(curpos+nbdofs); + for (int i=0;i<normaldofs;i++) vals.push_back(valsd[i]); + } +} + +template <class T> void NodeEnrichedFS<T>::f(MElement *ele, double u, double v, double w, std::vector<ValType> &vals) +{ + // We need parent parameters + MElement * elep; + if (ele->getParent()) elep = ele->getParent(); + else elep = ele; + + std::vector<ValType> valsd; + _SupportFS->f(elep,u,v,w,valsd); + int normaldofs=valsd.size(); + + // nbdofs element calcul + int nbdofs = normaldofs; + std::vector<int> EnrichedVertex; + for (int i=0 ;i<elep->getNumVertices();i++) + { + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(elep->getVertex(i)->getNum()); + if (it!=_TagEnrichedVertex->end()) + { + EnrichedVertex.push_back(i); + nbdofs = nbdofs + 1*_EnrichComp->size(); // enriched dof + } + } + + int curpos=vals.size(); + vals.resize(curpos+nbdofs); + + // first normal dofs + for (int i=0;i<normaldofs;i++) vals.push_back(valsd[i]); + + // then enriched dofs so the order is ex:(u1x,u2x,u3x,u1y,u2y,u3y,a2x,a2y,a3x,a3y) + if (nbdofs>normaldofs) // if enriched + { + // Enrichment function calcul + SPoint3 p; + elep->pnt(u, v, w, p); + double func; + func = (*_funcEnrichment)(p.x(), p.y(), p.z()); + + for (int i=0 ;i<EnrichedVertex.size();i++) + { + for (int j=0;j<_EnrichComp->size();j++) vals.push_back(valsd[EnrichedVertex[i]+(_EnrichComp->at(j))*elep->getNumVertices()]*func); + } + } +} + +template <class T> void NodeEnrichedFS<T>::gradf(MElement *ele, double u, double v, double w,std::vector<GradType> &grads) +{ + // We need parent parameters + MElement * elep; + if (ele->getParent()) elep = ele->getParent(); + else elep = ele; + + std::vector<GradType> gradsd; + _SupportFS->gradf(elep,u,v,w,gradsd); + + int normaldofs=gradsd.size(); + int nbdofs = normaldofs; + + GradType val; + + std::vector<int> EnrichedVertex; + + for (int i=0 ;i<elep->getNumVertices();i++) + { + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(elep->getVertex(i)->getNum()); + if (it!=_TagEnrichedVertex->end()) + { + EnrichedVertex.push_back(i); + nbdofs = nbdofs + 1*_EnrichComp->size(); // enriched dof + } + } + + + int curpos=grads.size(); + grads.reserve(curpos+nbdofs); + + // first normal dofs + for (int i=0;i<normaldofs;i++) grads.push_back(gradsd[i]); + + // then enriched dofs so the order is ex:(u1x,u2x,u3x,u1y,u2y,u3y,a2x,a2y,a3x,a3y) + if (nbdofs>normaldofs) // if enriched + { + // Enrichment function calcul + SPoint3 p; + elep->pnt(u, v, w, p); + double func; + func = (*_funcEnrichment)(p.x(), p.y(), p.z()); + + for (int i=0 ;i<EnrichedVertex.size();i++) + { + for (int j=0;j<_EnrichComp->size();j++){ + grads.push_back(gradsd[EnrichedVertex[i]+(_EnrichComp->at(j))*elep->getNumVertices()]*func); + } + } + } +} + + +template <class T> int NodeEnrichedFS<T>::getNumKeys(MElement *ele) +{ + MElement * elep; + if (ele->getParent()) elep = ele->getParent(); + else elep = ele; + + int normaldofs = _SupportFS->getNumKeys(ele); + int nbdofs = normaldofs; + + for (int i=0 ;i<elep->getNumVertices();i++) + { + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(elep->getVertex(i)->getNum()); + if (it!=_TagEnrichedVertex->end()) + { + nbdofs = nbdofs + 1*_EnrichComp->size(); // enriched dof + } + } + return nbdofs; +} + +template <class T> void NodeEnrichedFS<T>::getKeys(MElement *ele, std::vector<Dof> &keys) +{ + MElement * elep; + if (ele->getParent()) elep = ele->getParent(); + else elep = ele; + + int normalk=_SupportFS->getNumKeys(elep); + + std::vector<Dof> bufk; + bufk.reserve(normalk); + _SupportFS->getKeys(elep,bufk); + int normaldofs=bufk.size(); + int nbdofs = normaldofs; + + int curpos=keys.size(); + + std::vector<int> EnrichedVertex; + + for (int i=0 ;i<elep->getNumVertices();i++) + { + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(elep->getVertex(i)->getNum()); + if (it!=_TagEnrichedVertex->end()) + { + EnrichedVertex.push_back(i); + nbdofs = nbdofs + 1*_EnrichComp->size(); // enriched dof + } + } + + keys.reserve(curpos+nbdofs); + + // get keys so the order is ex:(u1x,u2x,u3x,u1y,u2y,u3y,a2x,a2y,a3x,a3y) + + for (int i=0;i<bufk.size();i++) keys.push_back(bufk[i]); + + for (int i=0 ;i<EnrichedVertex.size();i++) + { + for (int j=0;j<_EnrichComp->size();j++) + { + int i1,i2; + Dof::getTwoIntsFromType(bufk[EnrichedVertex[i]+_EnrichComp->at(j)*elep->getNumVertices()].getType(), i1,i2); + keys.push_back(Dof(bufk[EnrichedVertex[i]+_EnrichComp->at(j)*elep->getNumVertices()].getEntity(),Dof::createTypeWithTwoInts(_EnrichComp->at(j),i1+1))); + } + } +} + +template <class T> int NodeEnrichedFS<T>::getNumKeys(MVertex *ver) +{ + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(ver->getNum()); + if (it!=_TagEnrichedVertex->end()) // if node enriched + { + return _SupportFS->getNumKeys(ver)+_EnrichComp->size(); + } + else + { + return _SupportFS->getNumKeys(ver); + } +} + + +template <class T> void NodeEnrichedFS<T>::getKeys(MVertex *ver, std::vector<Dof> &keys) +{ + int normalkeys= _SupportFS->getNumKeys(ver); + std::vector<Dof> bufk; + bufk.reserve(normalkeys); + _SupportFS->getKeys(ver,bufk); + int normaldofs=bufk.size(); + int nbdofs = normaldofs; + int curpos=keys.size(); + + + std::set<int>::iterator it; + it = _TagEnrichedVertex->find(ver->getNum()); + if (it!=_TagEnrichedVertex->end()) + nbdofs = nbdofs+_EnrichComp->size(); // enriched dof + + keys.reserve(curpos+nbdofs); + + // normal dofs + for (int i=0;i<bufk.size();i++) keys.push_back(bufk[i]); + + // enriched dofs + if (nbdofs>normaldofs) + { + for (int j=0;j<_EnrichComp->size();j++) + { + int i1,i2; + Dof::getTwoIntsFromType(bufk[_EnrichComp->at(j)].getType(), i1,i2); + keys.push_back(Dof(bufk[_EnrichComp->at(j)].getEntity(),Dof::createTypeWithTwoInts(_EnrichComp->at(j),i1+1))); + } + } +}; diff --git a/contrib/arc/Classes/NodeEnrichedFS.h b/contrib/arc/Classes/NodeEnrichedFS.h new file mode 100644 index 0000000000000000000000000000000000000000..5767e94ae7d9e6f707340d34381ede051d9ca149 --- /dev/null +++ b/contrib/arc/Classes/NodeEnrichedFS.h @@ -0,0 +1,74 @@ +// +// Description : From function space to node enriched composant function space +// +// +// Author: <Eric Bechet>::<Boris Sedji>, 01/2010 +// +// Copyright: See COPYING file that comes with this distribution +// +// + +#ifndef _NODEENRICHEDFS_H_ +#define _NODEENRICHEDFS_H_ + +#include "SPoint3.h" +#include "SVector3.h" +#include "STensor3.h" +#include <vector> +#include <iterator> +#include <iostream> +#include "Numeric.h" +#include "MElement.h" +#include "MVertex.h" +#include "dofManager.h" +#include "functionSpace.h" +#include "simpleFunction.h" +#include <set> + + +template<class T> +class NodeEnrichedFS : public FunctionSpace<T> +{ + +public: + + typedef typename TensorialTraits<T>::ValType ValType; + typedef typename TensorialTraits<T>::GradType GradType; + typedef typename TensorialTraits<T>::HessType HessType; + typedef typename TensorialTraits<T>::DivType DivType; + typedef typename TensorialTraits<T>::CurlType CurlType; + +protected: + + std::set<int > *_TagEnrichedVertex; + std::vector<int> * _EnrichComp; + simpleFunction<double> *_funcEnrichment; + FunctionSpace<T> *_SupportFS; + +public: + + NodeEnrichedFS(FunctionSpace<T> * SupportFS, std::set<int > * TagEnrichedVertex , std::vector<int> * EnrichComp, simpleFunction<double> *funcEnrichment) + { + _SupportFS = SupportFS; + _TagEnrichedVertex = TagEnrichedVertex; + _funcEnrichment = funcEnrichment; + _EnrichComp = EnrichComp; + } + virtual ~NodeEnrichedFS() {} + // Shape function value in element at uvw + virtual void f(MElement *ele, double u, double v, double w, std::vector<ValType> &vals); + // Shape function value at vertex + virtual void f(MVertex *ver, std::vector<ValType> &vals) ; + // Grad Shape function value in element at uvw + virtual void gradf(MElement *ele, double u, double v, double w,std::vector<GradType> &grads) ; + // Shape function number for element + virtual int getNumKeys(MElement *ele); + // Dof keys for the element + virtual void getKeys(MElement *ele, std::vector<Dof> &keys) ; + // Shape function Number associate with the vertex + virtual int getNumKeys(MVertex *ver); + // Get Dof keys for the vertex + virtual void getKeys(MVertex *ver, std::vector<Dof> &keys); +}; + +#endif diff --git a/contrib/arc/Classes/VectorXFEMFS.cpp b/contrib/arc/Classes/VectorXFEMFS.cpp index 937c76ee4c8618160c1f14efc5fc645ba4e39b6e..268f8f26d1a5f7ef67e86a2b644aa23f83de2884 100644 --- a/contrib/arc/Classes/VectorXFEMFS.cpp +++ b/contrib/arc/Classes/VectorXFEMFS.cpp @@ -85,7 +85,10 @@ void ScalarLagrangeToXfemFS::gradf(MElement *ele, double u, double v, double w,s { std::set<int>::iterator it; it = _TagEnrichedVertex->find(ele->getVertex(i)->getNum()); - if (it!=_TagEnrichedVertex->end()){ndofs = ndofs + 1;} // enriched dof + if (it!=_TagEnrichedVertex->end()) + { + ndofs = ndofs + 1; // enriched dof + } } int curpos = grads.size(); diff --git a/contrib/arc/Classes/XFEMelasticitySolver.cpp b/contrib/arc/Classes/XFEMelasticitySolver.cpp index 5a45cdb0d8fee832e833fb955612bdf8056063ad..63e94e0b08617883e207b3b87980c13ae3cded98 100644 --- a/contrib/arc/Classes/XFEMelasticitySolver.cpp +++ b/contrib/arc/Classes/XFEMelasticitySolver.cpp @@ -34,6 +34,7 @@ #include "XFEMelasticitySolver.h" #include "FuncHeaviside.h" +#include "NodeEnrichedFS.cpp" XFEMelasticitySolver::~XFEMelasticitySolver() { @@ -78,16 +79,27 @@ void XFEMelasticitySolver::solve(){ } } + // enriched composant + _EnrichComp.push_back(0); + _EnrichComp.push_back(1); + //_EnrichComp.push_back(2); + + // level set definition (in .dat ??) double a(0.), b(1.), c(0.), d(-4.7); int n(1); // pour level set gLevelset *ls = new gLevelsetPlane(a, b, c, d, n); _funcEnrichment = new FuncHeaviside(ls); + // space function definition + FunctionSpace<SVector3> *NLagSpace; if (LagSpace) delete LagSpace; - if (_dim==3) LagSpace=new ScalarXFEMToVectorFS(_tag,_TagEnrichedVertex,_funcEnrichment); - if (_dim==2) LagSpace=new ScalarXFEMToVectorFS(_tag,ScalarXFEMToVectorFS::VECTOR_X,ScalarXFEMToVectorFS::VECTOR_Y,_TagEnrichedVertex,_funcEnrichment); + if (_dim==3) NLagSpace=new VectorLagrangeFunctionSpace(_tag); + if (_dim==2) NLagSpace=new VectorLagrangeFunctionSpace(_tag,VectorLagrangeFunctionSpace::VECTOR_X,VectorLagrangeFunctionSpace::VECTOR_Y); + + LagSpace = new NodeEnrichedFS<SVector3>(NLagSpace, &_TagEnrichedVertex ,&_EnrichComp, _funcEnrichment); + //LagSpace = new ScalarXFEMToVectorFS(_tag,ScalarXFEMToVectorFS::VECTOR_X,ScalarXFEMToVectorFS::VECTOR_Y,_TagEnrichedVertex,_funcEnrichment); // we first do all fixations. the behavior of the dofManager is to // give priority to fixations : when a dof is fixed, it cannot be diff --git a/contrib/arc/Classes/XFEMelasticitySolver.h b/contrib/arc/Classes/XFEMelasticitySolver.h index 408947f70f7cfcaf1ba0ec2eecac49bef07f0702..57266017f95091e9e0466d55c1035844b2df5822 100644 --- a/contrib/arc/Classes/XFEMelasticitySolver.h +++ b/contrib/arc/Classes/XFEMelasticitySolver.h @@ -19,6 +19,8 @@ class XFEMelasticitySolver : public elasticitySolver protected : // map containing the tag of vertex and enriched status std::set<int > _TagEnrichedVertex; + // enriched comp + std::vector<int> _EnrichComp; // simple multiplying function enrichment to enrich the space function simpleFunction <double> *_funcEnrichment;