From 1d0d72503ea0ec951a8bc86bd1c2d9193e752a28 Mon Sep 17 00:00:00 2001
From: Tuomas Karna <tuomas.karna@uclouvain.be>
Date: Fri, 4 Jun 2010 14:53:10 +0000
Subject: [PATCH] improved mesh extrusion in Lua: physical tags, partitioning,
second order prisms
---
Geo/GEntity.cpp | 9 ++++++++-
Geo/GEntity.h | 11 ++++++++++-
Geo/GModel.cpp | 10 ++++++++++
Geo/MElement.cpp | 3 +++
Geo/MLine.cpp | 13 ++++++++++---
Geo/MPrism.cpp | 16 ++++++++++++++++
Geo/MQuadrangle.cpp | 13 +++++++++++++
Geo/MTriangle.cpp | 16 ++++++++++++++++
Geo/MVertex.cpp | 5 +++++
Geo/MVertex.h | 1 +
Numeric/polynomialBasis.cpp | 27 +++++++++++++++++++++------
11 files changed, 113 insertions(+), 11 deletions(-)
diff --git a/Geo/GEntity.cpp b/Geo/GEntity.cpp
index 980f48f5b7..7ed9ad984f 100644
--- a/Geo/GEntity.cpp
+++ b/Geo/GEntity.cpp
@@ -102,5 +102,12 @@ void GEntity::registerBindings(binding *b)
mb->setDescription("do a dynamic cast of the GEntity to a GRegion (0 if wrong cast).");
mb = cb->addMethod("tag", &GEntity::tag);
mb->setDescription("return the tag of this entity.");
-
+ mb = cb->addMethod("getPhysicalEntities", &GEntity::getPhysicalEntities);
+ mb->setDescription("return a vector of all physical entities that this entity belongs to.");
+ mb = cb->addMethod("addPhysicalEntity", &GEntity::addPhysicalEntity);
+ mb->setArgNames("physicalGroupId",NULL);
+ mb->setDescription("add this element to a physical group.");
}
+
+
+
diff --git a/Geo/GEntity.h b/Geo/GEntity.h
index f446cf6322..9887aa4a7a 100644
--- a/Geo/GEntity.h
+++ b/Geo/GEntity.h
@@ -221,6 +221,15 @@ class GEntity {
int tag() const { return _tag; }
void setTag(int tag) { _tag = tag; }
+ // get/set physical entities
+ virtual void addPhysicalEntity(int physicalTag) {
+ physicals.push_back(physicalTag);
+ }
+ virtual std::vector<int> getPhysicalEntities() {
+ return physicals;
+ }
+
+
// returns the tag of the entity that its master entity (for mesh)
int meshMaster() const { return _meshMaster; }
void setMeshMaster(int m) { _meshMaster = m; }
@@ -234,7 +243,7 @@ class GEntity {
// get/set the visibility flag
virtual char getVisibility();
virtual void setVisibility(char val, bool recursive=false){ _visible = val; }
-
+
// get/set the selection flag
virtual char getSelection(){ return _selection; }
virtual void setSelection(char val){ _selection = val; }
diff --git a/Geo/GModel.cpp b/Geo/GModel.cpp
index 6e50e18961..9ca315801b 100644
--- a/Geo/GModel.cpp
+++ b/Geo/GModel.cpp
@@ -2131,4 +2131,14 @@ void GModel::registerBindings(binding *b)
cm->setDescription("set the model as the current (active) one");
cm = cb->setConstructor<GModel>();
cm->setDescription("Create an empty GModel");
+
+ cm = cb->addMethod("getPhysicalName", &GModel::getPhysicalName);
+ cm->setDescription("get the name of an physical group, identified by its dimension and number. "
+ "Returns empty string if physical name is not assigned");
+ cm->setArgNames("dim","number",NULL);
+ cm = cb->addMethod("setPhysicalName", &GModel::setPhysicalName);
+ cm->setDescription("set the name of an physical group, identified by its dimension and number. "
+ "If number=0, the first free number is chosen. Returns the number.");
+ cm->setArgNames("physicalName","dim","number",NULL);
+
}
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 66fd267c4e..b62648a362 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -1023,6 +1023,9 @@ void MElement::registerBindings(binding *b)
cm->setDescription("get the gmsh type of the element");
cm = cb->addMethod("getPartition", &MElement::getPartition);
cm->setDescription("get the partition to which the element belongs");
+ cm = cb->addMethod("setPartition", &MElement::setPartition);
+ cm->setDescription("set the partition to which the element belongs");
+ cm->setArgNames("iPartition",NULL);
cm = cb->addMethod("getPolynomialOrder", &MElement::getPolynomialOrder);
cm->setDescription("return the polynomial order the element");
cm = cb->addMethod("getDim", &MElement::getDim);
diff --git a/Geo/MLine.cpp b/Geo/MLine.cpp
index a8fcf4d930..efdd93e4b3 100644
--- a/Geo/MLine.cpp
+++ b/Geo/MLine.cpp
@@ -41,6 +41,9 @@ double MLine::getInnerRadius()
}
#include "Bindings.h"
+static MLine* MLine_binding(std::vector<MVertex*> v) {
+ return new MLine(v);
+}
void MLine::registerBindings(binding *b)
{
@@ -48,8 +51,12 @@ void MLine::registerBindings(binding *b)
cb->setDescription("A line mesh element.");
methodBinding *cm;
- cm = cb->setConstructor<MLine,MVertex*,MVertex*>();
- cm->setArgNames("v0","v1", NULL);
- cm->setDescription("Create a new line mesh element between v0 and v1.");
+ cm = cb->addMethod("MLine",&MLine_binding);
+ cm->setArgNames("vectorOfVertices", NULL);
+ cm->setDescription("Create a new line mesh element with the given vertices. "
+ "First 2 vertices must correspond to the beginning/end of the line.");
+// cm = cb->setConstructor<MLine,MVertex*,MVertex*>();
+// cm->setArgNames("v0","v1", NULL);
+// cm->setDescription("Create a new line mesh element between v0 and v1.");
cb->setParentClass<MElement>();
}
diff --git a/Geo/MPrism.cpp b/Geo/MPrism.cpp
index 5627aeeb10..e245eaeefc 100644
--- a/Geo/MPrism.cpp
+++ b/Geo/MPrism.cpp
@@ -126,6 +126,9 @@ void MPrism::getFaceInfo(const MFace &face, int &ithFace, int &sign, int &rot) c
Msg::Error("Could not get face information for prism %d", getNum());
}
#include "Bindings.h"
+static MPrism18* MPrism18_binding(std::vector<MVertex*> v) {
+ return new MPrism18(v);
+}
void MPrism::registerBindings(binding *b)
{
@@ -135,6 +138,19 @@ void MPrism::registerBindings(binding *b)
cm = cb->setConstructor<MPrism,MVertex*,MVertex*,MVertex*,MVertex*, MVertex*, MVertex*>();
cm->setArgNames("v0", "v1", "v2", "v3","v4","v5", NULL);
cm->setDescription("Create a new prism with top triangle (v0,v1,v2) and bottom one (v3,v4,v5).");
+ cm = cb->addMethod("getVolumeSign",&MPrism::getVolumeSign);
+ cm->setDescription("computes the sign of the element volume");
+ cm = cb->addMethod("revert",&MPrism::revert);
+ cm->setDescription("reorganises the element vertices so that volume is positive");
+
cb->setParentClass<MElement>();
+
+ cb = b->addClass<MPrism18>("MPrism18");
+ cb->setDescription("A mesh second-order prism.");
+ cm = cb->addMethod("MPrism18",&MPrism18_binding);
+// cm = cb->setConstructor<MPrism18_binding,std::vector<MVertex*> >();
+ cm->setArgNames("vectorOfVertices", NULL);
+ cm->setDescription("Create a new prism with vertices in vectorV (length=18).");
+ cb->setParentClass<MPrism>();
}
diff --git a/Geo/MQuadrangle.cpp b/Geo/MQuadrangle.cpp
index 35dec5246e..e5d4afe433 100644
--- a/Geo/MQuadrangle.cpp
+++ b/Geo/MQuadrangle.cpp
@@ -274,6 +274,11 @@ double MQuadrangle::getInnerRadius()
}
#include "Bindings.h"
+static MQuadrangle9* MQuadrangle9_binding(std::vector<MVertex*> v) {
+ return new MQuadrangle9(v);
+}
+
+
void MQuadrangle::registerBindings(binding *b)
{
classBinding *cb = b->addClass<MQuadrangle>("MQuadrangle");
@@ -283,4 +288,12 @@ void MQuadrangle::registerBindings(binding *b)
cm->setArgNames("v0", "v1", "v2", "v3", NULL);
cm->setDescription("Create a new quadrangle with vertices (v0,v1,v2,v3).");
cb->setParentClass<MElement>();
+
+ cb = b->addClass<MQuadrangle9>("MQuadrangle9");
+ cb->setDescription("A mesh second-order quadrangle.");
+ cm = cb->addMethod("MQuadrangle9",&MQuadrangle9_binding);
+// cm = cb->setConstructor<MQuadrangle9_binding,std::vector<MVertex*> >();
+ cm->setArgNames("vectorOfVertices", NULL);
+ cm->setDescription("Create a new quadrangle with vertices in vectorV (length=9).");
+ cb->setParentClass<MQuadrangle>();
}
diff --git a/Geo/MTriangle.cpp b/Geo/MTriangle.cpp
index fc65008a33..ce58c06cec 100644
--- a/Geo/MTriangle.cpp
+++ b/Geo/MTriangle.cpp
@@ -233,6 +233,9 @@ void MTriangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts)
*pts = getGQTPts(pOrder);
}
#include "Bindings.h"
+static MTriangle6* MTriangle6_binding(std::vector<MVertex*> v) {
+ return new MTriangle6(v);
+}
void MTriangle::registerBindings(binding *b)
{
@@ -243,5 +246,18 @@ void MTriangle::registerBindings(binding *b)
cm->setArgNames("v0", "v1", "v2", NULL);
cm->setDescription("Create a new triangle with vertices (v0,v1,v2).");
cb->setParentClass<MElement>();
+
+ cb = b->addClass<MTriangle6>("MTriangle6");
+ cb->setDescription("A mesh second-order triangle.");
+ cm = cb->addMethod("MTriangle6",&MTriangle6_binding);
+ cm->setArgNames("vectorOfVertices", NULL);
+ cm->setDescription("Create a new triangle with vertices given in the vector (length = 6).");
+ cb->setParentClass<MTriangle>();
+
+/* cb->setDescription("A mesh second-order triangle.");
+ cm = cb->setConstructor<MTriangle6_binding,std::vector<MVertex*> >();
+ cm->setArgNames("vectorOfVertices", NULL);
+ cm->setDescription("Create a new triangle with vertices given in the vector (length = 6).");
+ cb->setParentClass<MTriangle>();*/
}
diff --git a/Geo/MVertex.cpp b/Geo/MVertex.cpp
index bf523e1967..314c69861a 100644
--- a/Geo/MVertex.cpp
+++ b/Geo/MVertex.cpp
@@ -422,4 +422,9 @@ void MVertex::registerBindings(binding *b)
cm->setDescription("Create a new mesh vertex at (x,y,z).");
cm = cb->addMethod("getNum", &MVertex::getNum);
cm->setDescription("return the invariant vertex id");
+ cm = cb->addMethod("getPolynomialOrder", &MVertex::getPolynomialOrder);
+ cm->setDescription("return the polynomial order of vertex");
+ cm = cb->addMethod("setPolynomialOrder", &MVertex::setPolynomialOrder_binding);
+ cm->setDescription("assign the polynomial order of vertex");
+ cm->setArgNames("order",NULL);
}
diff --git a/Geo/MVertex.h b/Geo/MVertex.h
index b817b8bb68..12abb066e4 100644
--- a/Geo/MVertex.h
+++ b/Geo/MVertex.h
@@ -59,6 +59,7 @@ class MVertex{
// get the "polynomial order" of the vertex
inline int getPolynomialOrder(){ return _order; }
inline void setPolynomialOrder(char order){ _order = order; }
+ inline void setPolynomialOrder_binding(int order){ _order = order; }
// get/set the coordinates
inline double x() const { return _x; }
diff --git a/Numeric/polynomialBasis.cpp b/Numeric/polynomialBasis.cpp
index 39b93f36d2..79a94761f0 100644
--- a/Numeric/polynomialBasis.cpp
+++ b/Numeric/polynomialBasis.cpp
@@ -240,17 +240,32 @@ static fullMatrix<double> generatePascalPrism(int order)
int nbMonomials = (order + 1) * (order + 1) * (order + 2) / 2;
fullMatrix<double> monomials(nbMonomials, 3);
-
int index = 0;
fullMatrix<double> lineMonoms = generate1DMonomials(order);
fullMatrix<double> triMonoms = generatePascalTriangle(order);
- for (int j = 0; j < lineMonoms.size1(); j++) {
- for (int i = 0; i < triMonoms.size1(); i++) {
- monomials(index,0) = triMonoms(i,0);
- monomials(index,1) = triMonoms(i,1);
- monomials(index,2) = lineMonoms(j,0);
+ // store monomials in right order
+ for (int currentOrder = 0; currentOrder <= order; currentOrder++) {
+ int orderT = currentOrder, orderL = currentOrder;
+ for (orderL = 0; orderL < currentOrder; orderL++) {
+ // do all permutations of monoms for orderL, orderT
+ int iL = orderL;
+ for (int iT = (orderT)*(orderT+1)/2; iT < (orderT+1)*(orderT+2)/2 ;iT++) {
+ monomials(index,0) = triMonoms(iT,0);
+ monomials(index,1) = triMonoms(iT,1);
+ monomials(index,2) = lineMonoms(iL,0);
index ++;
+ }
}
+ orderL = currentOrder;
+ for (orderT = 0; orderT <= currentOrder; orderT++) {
+ int iL = orderL;
+ for (int iT = (orderT)*(orderT+1)/2; iT < (orderT+1)*(orderT+2)/2 ;iT++) {
+ monomials(index,0) = triMonoms(iT,0);
+ monomials(index,1) = triMonoms(iT,1);
+ monomials(index,2) = lineMonoms(iL,0);
+ index ++;
+ }
+ }
}
// monomials.print("Pri monoms");
return monomials;
--
GitLab