From 4a8921d7c8f5db2ead52902826c9da13b54f16c9 Mon Sep 17 00:00:00 2001
From: Amaury Johnen <amaury.johnen@uclouvain.be>
Date: Fri, 11 Aug 2017 12:24:57 +0200
Subject: [PATCH] give interpolation matrices
---
Plugin/AnalyseCurvedMesh.cpp | 75 ++++++++++++++++++++++--------------
Plugin/AnalyseCurvedMesh.h | 3 ++
Post/PViewData.cpp | 5 +++
Post/PViewData.h | 1 +
4 files changed, 55 insertions(+), 29 deletions(-)
diff --git a/Plugin/AnalyseCurvedMesh.cpp b/Plugin/AnalyseCurvedMesh.cpp
index 8343055c84..a0becda9c6 100644
--- a/Plugin/AnalyseCurvedMesh.cpp
+++ b/Plugin/AnalyseCurvedMesh.cpp
@@ -25,6 +25,7 @@
#include "MHexahedron.h"
#include "MTetrahedron.h"
#include "MPrism.h"
+#include "BasisFactory.h"
class bezierBasis;
@@ -34,10 +35,10 @@ StringXNumber CurvedMeshOptions_Number[] = {
{GMSH_FULLRC, "ICN measure", NULL, 1},
{GMSH_FULLRC, "Hiding threshold", NULL, 9},
{GMSH_FULLRC, "Draw PView", NULL, 0},
- {GMSH_FULLRC, "Recompute", NULL, 0},
+ {GMSH_FULLRC, "Recompute", NULL, 1},
{GMSH_FULLRC, "Dimension of elements", NULL, -1},
//{GMSH_FULLRC, "Element to print quality", NULL, 2485}
- {GMSH_FULLRC, "Element to print quality", NULL, 2555}
+ {GMSH_FULLRC, "Element to print quality", NULL, 2901}
};
extern "C"
@@ -181,7 +182,15 @@ PView* GMSH_AnalyseCurvedMeshPlugin::execute(PView *v)
dataPVelement[_hoElement->getNum()] = _jacElementToScan;
std::stringstream name;
name << "Jacobian elem " << _numElementToScan;
- new PView(name.str().c_str(), "ElementNodeData", _m, dataPVelement, 0, 1);
+ PView *view = new PView(name.str().c_str(), "ElementNodeData",
+ _m, dataPVelement, 0, 1);
+ const nodalBasis *fs = BasisFactory::getNodalBasis(_hoElement->getTypeForMSH());
+ const polynomialBasis *pfs = dynamic_cast<const polynomialBasis*>(fs);
+ PViewData *viewData = view->getData();
+ viewData->deleteInterpolationMatrices(_hoElement->getType());
+ viewData->setInterpolationMatrices(_hoElement->getType(),
+ pfs->coefficients, pfs->monomials,
+ pfs->coefficients, pfs->monomials);
}
// Create PView
@@ -364,32 +373,7 @@ void GMSH_AnalyseCurvedMeshPlugin::_computeMinMaxJandValidity(int dim)
if (min < 0 && max < 0) ++cntInverted;
progress.next();
- if (el->getNum() == _numElementToScan) {
- _elementToScan = el;
- fullMatrix<double> points =
- _elementToScan->getFunctionSpace(_viewOrder)->points;
- int tag = ElementType::getTag(_elementToScan->getType(), _viewOrder);
- std::vector<MVertex *> v;
- for (int k = 0; k < points.size1(); k++) {
- double t1 = points(k, 0);
- double t2 = points(k, 1);
- double t3 = points(k, 2);
- SPoint3 pos;
- _elementToScan->pnt(t1, t2, t3, pos);
- v.push_back(new MVertex(pos.x(), pos.y(), pos.z()));
- }
- MElementFactory factory;
- _hoElement = factory.create(tag, v);
-
- _addElementInEntity(_hoElement, entity);
-
- fullVector<double> jac;
- jacobianBasedQuality::sampleJacobian(_elementToScan, _viewOrder,
- jac, normals);
- for (int j = 0; j < jac.size(); ++j) {
- _jacElementToScan.push_back(jac(j));
- }
- }
+ _computeJacobianToScan(el, normals);
}
delete normals;
}
@@ -563,6 +547,39 @@ void GMSH_AnalyseCurvedMeshPlugin::_printStatICN()
infminI, avgminI, supminI);
}
+void GMSH_AnalyseCurvedMeshPlugin::_computeJacobianToScan(MElement *el,
+ GEntity *entity,
+ const fullMatrix<double> *normals)
+{
+ if (el->getNum() != _numElementToScan) return;
+
+ _elementToScan = el;
+ fullMatrix<double> points =
+ _elementToScan->getFunctionSpace(_viewOrder)->points;
+ int tag = ElementType::getTag(_elementToScan->getType(), _viewOrder);
+ std::vector<MVertex *> v;
+ for (int k = 0; k < points.size1(); k++) {
+ double t1 = points(k, 0);
+ double t2 = points(k, 1);
+ double t3 = points(k, 2);
+ SPoint3 pos;
+ _elementToScan->pnt(t1, t2, t3, pos);
+ v.push_back(new MVertex(pos.x(), pos.y(), pos.z()));
+ }
+ MElementFactory factory;
+ _hoElement = factory.create(tag, v);
+
+ _addElementInEntity(_hoElement, entity);
+
+ fullVector<double> jac;
+ jacobianBasedQuality::sampleJacobian(_elementToScan, _viewOrder,
+ jac, normals);
+ _jacElementToScan.clear();
+ for (int j = 0; j < jac.size(); ++j) {
+ _jacElementToScan.push_back(jac(j));
+ }
+}
+
void GMSH_AnalyseCurvedMeshPlugin::_addElementInEntity(MElement *element,
GEntity *entity)
{
diff --git a/Plugin/AnalyseCurvedMesh.h b/Plugin/AnalyseCurvedMesh.h
index 7722695b87..51d572a114 100644
--- a/Plugin/AnalyseCurvedMesh.h
+++ b/Plugin/AnalyseCurvedMesh.h
@@ -87,6 +87,9 @@ private :
void _printStatJacobian();
void _printStatIGE();
void _printStatICN();
+
+ void _computeJacobianToScan(MElement*, GEntity*,
+ const fullMatrix<double> *normals);
void _addElementInEntity(MElement*, GEntity*);
};
diff --git a/Post/PViewData.cpp b/Post/PViewData.cpp
index 59d56e5155..0473d24ec2 100644
--- a/Post/PViewData.cpp
+++ b/Post/PViewData.cpp
@@ -184,6 +184,11 @@ bool PViewData::haveInterpolationMatrices(int type)
return _interpolation.count(type) ? true : false;
}
+void PViewData::deleteInterpolationMatrices(int type)
+{
+ _interpolation.erase(type);
+}
+
void PViewData::removeInterpolationScheme(const std::string &name)
{
std::map<std::string, interpolationMatrices>::iterator it = _interpolationSchemes.find(name);
diff --git a/Post/PViewData.h b/Post/PViewData.h
index 9b5449ec12..922083f180 100644
--- a/Post/PViewData.h
+++ b/Post/PViewData.h
@@ -225,6 +225,7 @@ class PViewData {
const fullMatrix<double> &expGeo);
int getInterpolationMatrices(int type, std::vector<fullMatrix<double>*> &p);
bool haveInterpolationMatrices(int type=0);
+ void deleteInterpolationMatrices(int type=0);
// access to global interpolation schemes
static void removeInterpolationScheme(const std::string &name);
--
GitLab