diff --git a/FunctionSpace/CMakeLists.txt b/FunctionSpace/CMakeLists.txt index 6ba65dfec500f1b32fd8e20066f60a236e8c4432..724e30c5ce17c4ede7b338f26005e8811178095c 100644 --- a/FunctionSpace/CMakeLists.txt +++ b/FunctionSpace/CMakeLists.txt @@ -14,7 +14,7 @@ set(SRC TriReferenceSpace.cpp QuadReferenceSpace.cpp TetReferenceSpace.cpp - #HexReferenceSpace.cpp + HexReferenceSpace.cpp ReferenceSpaceLagrange.cpp TriLagrangeReferenceSpace.cpp diff --git a/FunctionSpace/HexReferenceSpace.cpp b/FunctionSpace/HexReferenceSpace.cpp new file mode 100644 index 0000000000000000000000000000000000000000..ce427cdd6076b7e3e05a1039e4a79528104e22c9 --- /dev/null +++ b/FunctionSpace/HexReferenceSpace.cpp @@ -0,0 +1,89 @@ +#include <sstream> +#include "HexReferenceSpace.h" +#include "MHexahedron.h" + +using namespace std; + +HexReferenceSpace::HexReferenceSpace(void){ + // Vertex Definition // + nVertex = 8; + + // Edge Definition // + const size_t nEdge = 12; + refEdgeNodeIdx.resize(nEdge); + + for(size_t i = 0; i < nEdge; i++){ + refEdgeNodeIdx[i].resize(2); // Two Nodes per Edge + refEdgeNodeIdx[i][0] = MHexahedron::edges_hexa(i, 0); + refEdgeNodeIdx[i][1] = MHexahedron::edges_hexa(i, 1); + } + + // Face Definition // + size_t nFace = 6; + refFaceNodeIdx.resize(nFace); + + for(size_t i = 0; i < nFace; i++){ + refFaceNodeIdx[i].resize(4); // Four Nodes per Face + refFaceNodeIdx[i][0] = MHexahedron::faces_hexa(i, 0); + refFaceNodeIdx[i][1] = MHexahedron::faces_hexa(i, 1); + refFaceNodeIdx[i][2] = MHexahedron::faces_hexa(i, 2); + refFaceNodeIdx[i][3] = MHexahedron::faces_hexa(i, 3); + } + + // Init All // + init(); +} + +HexReferenceSpace::~HexReferenceSpace(void){ +} + +string HexReferenceSpace::toLatex(void) const{ + //const size_t nRefSpace = refSpaceNodeId.size(); + stringstream stream; + /* + stream << "\\documentclass{article}" << endl << endl + + << "\\usepackage{longtable}" << endl + << "\\usepackage{tikz}" << endl + << "\\usetikzlibrary{arrows}" << endl << endl + + << "\\begin{document}" << endl + << "\\tikzstyle{vertex} = [circle, fill = black!25]" << endl + << "\\tikzstyle{line} = [draw, thick, black, -latex']" << endl + << endl + + << "\\begin{longtable}{ccc}" << endl << endl; + + for(size_t s = 0; s < nRefSpace; s++){ + stream << "\\begin{tikzpicture}" << endl + + << "\\node[vertex] (n0) at(0, 0) {$" << refSpaceNodeId[s][0] << "$};" + << endl + << "\\node[vertex] (n1) at(3, 0) {$" << refSpaceNodeId[s][1] << "$};" + << endl + << "\\node[vertex] (n2) at(0, 3) {$" << refSpaceNodeId[s][2] << "$};" + << endl + << "\\node[vertex] (n3) at(1, 1) {$" << refSpaceNodeId[s][3] << "$};" + << endl + << endl; + + for(size_t e = 0; e < 6; e++) + stream << "\\path[line]" + << " (n" << orderedEdgeNodeIdx[s][e][0] << ")" + << " -- " + << " (n" << orderedEdgeNodeIdx[s][e][1] << ");" + << endl; + + if((s + 1) % 3) + stream << "\\end{tikzpicture} & " << endl << endl; + + else + stream << "\\end{tikzpicture} \\\\ \\\\" << endl << endl; + } + + stream << "\\end{longtable}" << endl + << "\\end{document}" << endl; + */ + stream << "Not Implemented" << endl; + return stream.str(); +} diff --git a/FunctionSpace/HexReferenceSpace.h b/FunctionSpace/HexReferenceSpace.h new file mode 100644 index 0000000000000000000000000000000000000000..bd03576f047eaa6c3092806564c7fd424cc7c406 --- /dev/null +++ b/FunctionSpace/HexReferenceSpace.h @@ -0,0 +1,31 @@ +#ifndef _HEXREFERENCESPACE_H_ +#define _HEXREFERENCESPACE_H_ + +#include <string> +#include "ReferenceSpace.h" + +/** + @class HexReferenceSpace + @brief ReferenceSpace for Hexahedron + + This class implements a ReferenceSpace for a Hexahedron. + */ + +class HexReferenceSpace: public ReferenceSpace{ + public: + HexReferenceSpace(void); + virtual ~HexReferenceSpace(void); + + virtual std::string toLatex(void) const; +}; + +/** + @fn HexReferenceSpace::HexReferenceSpace + Instatiate a new ReferenceSpace for a Hexahedron + ** + + @fn HexReferenceSpace::~HexReferenceSpace + Deletes this HexReferenceSpace +*/ + +#endif diff --git a/FunctionSpace/LineReferenceSpace.cpp b/FunctionSpace/LineReferenceSpace.cpp index c87a3d56bdc2a4ca6568e1fe73e837db22308e29..7106c8b6f7230ffd29343d835b41ca598194372d 100644 --- a/FunctionSpace/LineReferenceSpace.cpp +++ b/FunctionSpace/LineReferenceSpace.cpp @@ -16,7 +16,6 @@ LineReferenceSpace::LineReferenceSpace(void){ // Face Definition // refFaceNodeIdx.clear(); // No Face in Line - parallelFaceId.clear(); // And no parallel Face init(); } diff --git a/FunctionSpace/PermutationTree.cpp b/FunctionSpace/PermutationTree.cpp index 7f26080a1817a458b68cb5b3736e4a393a6e0081..018ba7dd8b08381396dc0fe0cb4c5f416ab28961 100644 --- a/FunctionSpace/PermutationTree.cpp +++ b/FunctionSpace/PermutationTree.cpp @@ -1,3 +1,4 @@ +#include <map> #include <sstream> #include "Exception.h" @@ -141,6 +142,29 @@ void PermutationTree::fillWithPermutation(size_t permutationId, } } +std::vector<std::pair<size_t, size_t> > +PermutationTree::getAllTagsCount(void) const{ + // Init Temp map + map<size_t, size_t> counter; + + // Loop on leafs and populate map + const size_t nLeaf = leaf.size(); + pair<map<size_t, size_t>::iterator, bool> inserted; + + for(size_t i = 0; i < nLeaf; i++){ + // Try insert new tag with counter to one + inserted = counter.insert(std::pair<size_t, size_t>(leaf[i]->tag, 1)); + + if(!inserted.second) + // If known tag -- add one to mapped value + inserted.first->second++; + } + + // Serialize this map and return + std::vector<std::pair<size_t, size_t> > ret(counter.begin(), counter.end()); + return ret; +} + string PermutationTree::toString(void) const{ // Number of Permutation // const size_t max_t = -1; diff --git a/FunctionSpace/PermutationTree.h b/FunctionSpace/PermutationTree.h index 9c54628692ff50c819c85e19a808502b3ce5dd16..f0cd7e99a6f90036afc5b679b01f007faced0714 100644 --- a/FunctionSpace/PermutationTree.h +++ b/FunctionSpace/PermutationTree.h @@ -58,7 +58,7 @@ class PermutationTree{ void addTagToPermutation(size_t permutationId, size_t tag); size_t getTagFromPermutation(size_t permutationId); - //std::vector<size_t> getAllTags(void) const; + std::vector<std::pair<size_t, size_t> > getAllTagsCount(void) const; std::string toString(void) const; @@ -120,6 +120,14 @@ class PermutationTree{ @return Returns the tag of the given sequence ** + @fn PermutationTree::getAllTagsCount + @return Returns a vector of pair such that: + @li The first entry is a tag + @li The second entry is the number of node having this tag + + The returned vector has en entry for each possible tag + ** + @fn PermutationTree::toString @return Returns a string describing this PermutationTree */ diff --git a/FunctionSpace/QuadReferenceSpace.cpp b/FunctionSpace/QuadReferenceSpace.cpp index 318c3e3c9e2ac8bd7f33c1d1f0fe56c500b917ac..47495dd8d718e063b1e65e4d0f3b0d504ab6cf53 100644 --- a/FunctionSpace/QuadReferenceSpace.cpp +++ b/FunctionSpace/QuadReferenceSpace.cpp @@ -27,9 +27,6 @@ QuadReferenceSpace::QuadReferenceSpace(void){ refFaceNodeIdx[0][2] = 2; refFaceNodeIdx[0][3] = 3; - parallelFaceId.resize(1); - parallelFaceId[0] = -1; // Face zero has no parallel face - // Init All // init(); } diff --git a/FunctionSpace/ReferenceSpace.cpp b/FunctionSpace/ReferenceSpace.cpp index 3959ee9f12a00e088ca44d2980cf76af195fc420..e3a3bd3062138ea49d8f3ed991376a4628e10995 100644 --- a/FunctionSpace/ReferenceSpace.cpp +++ b/FunctionSpace/ReferenceSpace.cpp @@ -1,5 +1,6 @@ #include <algorithm> #include <sstream> +#include <set> #include "Exception.h" #include "Numeric.h" @@ -111,6 +112,7 @@ findCyclicPermutation(list<size_t>& listOfTrueReferenceSpace, // --> this Permutation is a new Reference Space if(!match.first){ listOfTrueReferenceSpace.push_back(i); + listOfRefNodeIndexPermutation.push_back(unPermutedIndex); listOfReverseNodeIndexPermutation.push_back(unPermutedIndex); @@ -121,22 +123,24 @@ findCyclicPermutation(list<size_t>& listOfTrueReferenceSpace, else{ listOfRefNodeIndexPermutation.push_back(match.second); listOfReverseNodeIndexPermutation.push_back(match.third); + pTree->addTagToPermutation(i, *it); } } } ReferenceSpace::triplet ReferenceSpace:: -isCyclicPermutation(vector<size_t>& pTest, - vector<size_t>& pRef){ +isCyclicPermutation(const vector<size_t>& pTest, + const vector<size_t>& pRef) const{ + + // Triplet to return + triplet tri; // If no Face, we have no Cyclic Permutation if(!refFaceNodeIdx.size()){ - triplet tri = { - false, - vector<size_t>(0), - vector<size_t>(0) - }; + tri.first = false; + tri.second = vector<size_t>(0); + tri.third = vector<size_t>(0); return tri; } @@ -153,11 +157,9 @@ isCyclicPermutation(vector<size_t>& pTest, // If no corresponding face found, we have no Cyclic Permutation if(testFaceId == (size_t)(-1)){ - triplet tri = { - false, - vector<size_t>(0), - vector<size_t>(0) - }; + tri.first = false; + tri.second = vector<size_t>(0); + tri.third = vector<size_t>(0); return tri; } @@ -169,18 +171,29 @@ isCyclicPermutation(vector<size_t>& pTest, for(size_t i = 0; i < nNodeInFaceTest; i++) testNode[i] = pTest[refFaceNodeIdx[testFaceId][i]]; - // Return Triplet - triplet tri = { - isFacePermutation(refNode, testNode), - getRefIndexPermutation(pRef, pTest), - getReverseIndexPermutation(pRef, pTest) - }; + // Test if we have a face cyclic permutation + size_t isCyclic = isFacePermutation(refNode, testNode); + + // Test ifwe have the same connectivity + bool isSameConnectivity = isSameEdge(pTest, pRef); + + if(isCyclic && isSameConnectivity){ + tri.first = true; + tri.second = getRefIndexPermutation(pRef, pTest); + tri.third = getReverseIndexPermutation(pRef, pTest); + } + + else{ + tri.first = false; + tri.second = vector<size_t>(0); + tri.third = vector<size_t>(0); + } return tri; } -size_t ReferenceSpace::findCorrespondingFace(vector<size_t>& face, - vector<size_t>& node){ +size_t ReferenceSpace::findCorrespondingFace(const vector<size_t>& face, + const vector<size_t>& node) const{ // Init Stuff // const size_t nFace = refFaceNodeIdx.size(); const size_t faceSize = face.size(); @@ -216,29 +229,33 @@ size_t ReferenceSpace::findCorrespondingFace(vector<size_t>& face, return -1; } -bool ReferenceSpace::haveSameNode(vector<size_t> face0, - vector<size_t> face1){ +bool ReferenceSpace::haveSameNode(const vector<size_t>& face0, + const vector<size_t>& face1){ // Check Sizes const size_t size = face0.size(); if(size != face1.size()) return false; - // Sort both vectors - sort(face0.begin(), face0.end()); - sort(face1.begin(), face1.end()); + // Copy vector + vector<size_t> cFace0(face0); + vector<size_t> cFace1(face1); + + // Sort both copy vectors + sort(cFace0.begin(), cFace0.end()); + sort(cFace1.begin(), cFace1.end()); // Check if elements are the same and at same position bool match = true; for(size_t i = 0; i < size && match; i++) - match = (face0[i] == face1[i]); + match = (cFace0[i] == cFace1[i]); return match; } -bool ReferenceSpace::isFacePermutation(vector<size_t>& refNode, - vector<size_t>& testNode){ +bool ReferenceSpace::isFacePermutation(const vector<size_t>& refNode, + const vector<size_t>& testNode){ // Check Size const size_t size = refNode.size(); @@ -264,8 +281,64 @@ bool ReferenceSpace::isFacePermutation(vector<size_t>& refNode, return match; } -vector<size_t> ReferenceSpace::getRefIndexPermutation(vector<size_t>& ref, - vector<size_t>& test){ +bool ReferenceSpace::isSameEdge(const std::vector<size_t>& pTest, + const std::vector<size_t>& pRef) const{ + // Set of Reference Edges + const size_t nEdge = refEdgeNodeIdx.size(); + std::set<std::vector<size_t>, EdgeComparator> refEdge; + std::vector<size_t> tmp(2); + + for(size_t e = 0; e < nEdge; e++){ + tmp[0] = pRef[refEdgeNodeIdx[e][0]]; + tmp[1] = pRef[refEdgeNodeIdx[e][1]]; + refEdge.insert(tmp); + } + + // Compare Test Edges + bool match = true; + for(size_t e = 0; e < nEdge && match; e++){ + tmp[0] = pTest[refEdgeNodeIdx[e][0]]; + tmp[1] = pTest[refEdgeNodeIdx[e][1]]; + + match = (refEdge.find(tmp) != refEdge.end()); + } + + return match; +} + +bool ReferenceSpace::EdgeComparator::operator()(const std::vector<size_t>& a, + const std::vector<size_t>& b){ + size_t maxA, minA; + size_t maxB, minB; + + if(a[0] < a[1]){ + maxA = a[1]; + minA = a[0]; + } + + else{ + maxA = a[0]; + minA = a[1]; + } + + if(b[0] < b[1]){ + maxB = b[1]; + minB = b[0]; + } + + else{ + maxB = b[0]; + minB = b[1]; + } + + return + ((minA != minB) && (minA < minB)) || + ((minA == minB) && (maxA < maxB)); +} + +vector<size_t> +ReferenceSpace::getRefIndexPermutation(const vector<size_t>& ref, + const vector<size_t>& test) const{ const size_t size = ref.size(); vector<size_t> idxVec(ref.size()); size_t idx; @@ -282,8 +355,9 @@ vector<size_t> ReferenceSpace::getRefIndexPermutation(vector<size_t>& ref, return idxVec; } -vector<size_t> ReferenceSpace::getReverseIndexPermutation(vector<size_t>& ref, - vector<size_t>& test){ +vector<size_t> +ReferenceSpace::getReverseIndexPermutation(const vector<size_t>& ref, + const vector<size_t>& test) const{ const size_t size = ref.size(); vector<size_t> idxVec(ref.size()); size_t idx; @@ -361,7 +435,7 @@ void ReferenceSpace::getOrderedFace(void){ refSpaceNodeId[s], orderedFaceNodeIdx[s][f]); if(nNodeInFace == 4) - correctQuadFaceNodeIdx(orderedFaceNodeIdx[s][f]); + correctQuadFaceNodeIdx(f, orderedFaceNodeIdx[s][f]); } } } @@ -390,15 +464,28 @@ orderRefEntityForGivenRefSpace(vector<size_t>& refEntityNodeIdx, } void ReferenceSpace:: -correctQuadFaceNodeIdx(vector<size_t>& correctedQuadFaceNodeIdx){ - +correctQuadFaceNodeIdx(size_t faceId, + vector<size_t>& correctedQuadFaceNodeIdx){ // Get : // // correctedQuadFaceNodeIdx[2] // // is *opposite* to // // correctedQuadFaceNodeIdx[0] // - size_t opposite = (correctedQuadFaceNodeIdx[0] + 2) % 4; + // Get opposite Node // + // Find node correctedQuadFaceNodeIdx[0] in refFaceNodeId[faceId] + size_t refIdx = 0; + while((refIdx < 4) + && + (refFaceNodeIdx[faceId][refIdx] != correctedQuadFaceNodeIdx[0])){ + + refIdx++; + } + + // Get opposite Node + size_t opposite = refFaceNodeIdx[faceId][(refIdx + 2) % 4]; + + // If correction needed, do it // if(correctedQuadFaceNodeIdx[2] != opposite){ // Find opposite node index size_t tmp; @@ -693,11 +780,7 @@ string ReferenceSpace::toString(void) const{ const size_t nFace = refFaceNodeIdx.size(); stringstream stream; - // Permutation Tree // - stream << "Permutation Tree:" << endl - << pTree->toString() << endl; - - // Reference Space // + // Reference Space Node IDs // stream << "Reference Space Node IDs:" << endl; for(size_t s = 0; s < nRefSpace; s++){ @@ -709,6 +792,7 @@ string ReferenceSpace::toString(void) const{ stream << refSpaceNodeId[s][nVertex - 1] << "]" << endl; } + // Oriented Edges Node Index // stream << endl << "Ordered Edge Node Index:" << endl; for(size_t s = 0; s < nRefSpace; s++){ @@ -722,6 +806,7 @@ string ReferenceSpace::toString(void) const{ << endl; } + // Oriented Faces Node Index // stream << endl << "Ordered Face Node Index:" << endl; for(size_t s = 0; s < nRefSpace; s++){ @@ -740,9 +825,6 @@ string ReferenceSpace::toString(void) const{ } } - //stream << endl << pTree->getNPermutation() << endl; - stream << endl; - return stream.str(); } diff --git a/FunctionSpace/ReferenceSpace.h b/FunctionSpace/ReferenceSpace.h index 910eefb59189416e4b667dbe92413daa0c4295f5..ebaa7becccf11ff111881b327e183701a2fce415 100644 --- a/FunctionSpace/ReferenceSpace.h +++ b/FunctionSpace/ReferenceSpace.h @@ -19,6 +19,13 @@ */ class ReferenceSpace{ + private: + class EdgeComparator{ + public: + bool operator()(const std::vector<size_t>& a, + const std::vector<size_t>& b); + }; + private: typedef struct{ bool first; @@ -31,7 +38,6 @@ class ReferenceSpace{ size_t nVertex; std::vector<std::vector<size_t> > refEdgeNodeIdx; std::vector<std::vector<size_t> > refFaceNodeIdx; - std::vector<size_t> parallelFaceId; // Permutation Tree // PermutationTree* pTree; @@ -101,23 +107,31 @@ class ReferenceSpace{ std::list<std::vector<size_t> >& listOfRefNodeIndexPermutation, std::list<std::vector<size_t> >& listOfReverseNodeIndexPermutation); - triplet isCyclicPermutation(std::vector<size_t>& pTest, - std::vector<size_t>& pRef); + triplet isCyclicPermutation(const std::vector<size_t>& pTest, + const std::vector<size_t>& pRef) const; + + size_t findCorrespondingFace(const std::vector<size_t>& face, + const std::vector<size_t>& node) const; - size_t findCorrespondingFace(std::vector<size_t>& face, - std::vector<size_t>& node); + static bool haveSameNode(const std::vector<size_t>& face0, + const std::vector<size_t>& face1); - static bool haveSameNode(std::vector<size_t> face0, - std::vector<size_t> face1); + static bool isFacePermutation(const std::vector<size_t>& refNode, + const std::vector<size_t>& testNode); - static bool isFacePermutation(std::vector<size_t>& refNode, - std::vector<size_t>& testNode); + bool isSameEdge(const std::vector<size_t>& pTest, + const std::vector<size_t>& pRef) const; - std::vector<size_t> getRefIndexPermutation(std::vector<size_t>& ref, - std::vector<size_t>& test); + static bool edgeComparator(const std::vector<size_t>& a, + const std::vector<size_t>& b); - std::vector<size_t> getReverseIndexPermutation(std::vector<size_t>& ref, - std::vector<size_t>& test); + std::vector<size_t> + getRefIndexPermutation(const std::vector<size_t>& ref, + const std::vector<size_t>& test) const; + + std::vector<size_t> + getReverseIndexPermutation(const std::vector<size_t>& ref, + const std::vector<size_t>& test) const; size_t getPermutationIdx(const MElement& element) const; @@ -126,8 +140,8 @@ class ReferenceSpace{ std::vector<size_t>& refSpaceNodeId, std::vector<size_t>& orderedEntityNodeIdx); - static void - correctQuadFaceNodeIdx(std::vector<size_t>& correctedQuadFaceNodeIdx); + void correctQuadFaceNodeIdx(size_t faceId, + std::vector<size_t>& correctedQuadFaceNodeIdx); static bool sortPredicate(const std::pair<size_t, size_t>& a, const std::pair<size_t, size_t>& b); diff --git a/FunctionSpace/TetReferenceSpace.cpp b/FunctionSpace/TetReferenceSpace.cpp index 8f5fbe33e60756def9ad180c86267f01c72bb008..402d9104df9feed0af1cb70314c4e2a7c2e7c7df 100644 --- a/FunctionSpace/TetReferenceSpace.cpp +++ b/FunctionSpace/TetReferenceSpace.cpp @@ -29,11 +29,6 @@ TetReferenceSpace::TetReferenceSpace(void){ refFaceNodeIdx[i][2] = MTetrahedron::faces_tetra(i, 2); } - parallelFaceId.resize(nFace); - - for(size_t i = 0; i < nFace; i++) - parallelFaceId[i] = -1; // Face i has no parallel face - // Init All // init(); } diff --git a/FunctionSpace/TriReferenceSpace.cpp b/FunctionSpace/TriReferenceSpace.cpp index 0494373e00b8ea5fae9e0b06d3d4fefd96054c16..e7b9399f6233fca7e7fe6fd624c97650226e5414 100644 --- a/FunctionSpace/TriReferenceSpace.cpp +++ b/FunctionSpace/TriReferenceSpace.cpp @@ -26,9 +26,6 @@ TriReferenceSpace::TriReferenceSpace(void){ refFaceNodeIdx[0][1] = 1; refFaceNodeIdx[0][2] = 2; - parallelFaceId.resize(1); - parallelFaceId[0] = -1; // Face zero has no parallel face - // Init All // init(); }