From b6f9dbabf8bb13fba8049277ea1cf937aac453d3 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Thu, 27 Nov 2008 22:19:27 +0000
Subject: [PATCH] fix tetra edge lighting fir order >1
---
Geo/MElement.cpp | 109 ++++++++++++++++++------------------------
doc/texinfo/gmsh.texi | 10 +++-
2 files changed, 55 insertions(+), 64 deletions(-)
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 7cb9df82d6..5b094eae01 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -813,13 +813,11 @@ const gmshFunctionSpace* MTriangle::getFunctionSpace(int o) const
return 0;
}
-#define NUM_SUB_EDGES (CTX.mesh.num_sub_edges)
+int MTriangleN::getNumEdgesRep(){ return 3 * CTX.mesh.num_sub_edges; }
+int MTriangle6::getNumEdgesRep(){ return 3 * CTX.mesh.num_sub_edges; }
-int MTriangleN::getNumEdgesRep(){ return 3 * NUM_SUB_EDGES; }
-int MTriangle6::getNumEdgesRep(){ return 3 * NUM_SUB_EDGES; }
-
-static void _myGetEdgeRep(MTriangle *t, int num, double *x, double *y, double *z, SVector3 *n,
- int numSubEdges)
+static void _myGetEdgeRep(MTriangle *t, int num, double *x, double *y, double *z,
+ SVector3 *n, int numSubEdges)
{
n[0] = n[1] = n[2] = t->getFace(0).normal();
@@ -853,23 +851,26 @@ static void _myGetEdgeRep(MTriangle *t, int num, double *x, double *y, double *z
}
}
-void MTriangleN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n){
- _myGetEdgeRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+void MTriangleN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+{
+ _myGetEdgeRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
-void MTriangle6::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n){
- _myGetEdgeRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+
+void MTriangle6::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
+{
+ _myGetEdgeRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
-int MTriangle6::getNumFacesRep(){ return NUM_SUB_EDGES * NUM_SUB_EDGES; }
-int MTriangleN::getNumFacesRep(){ return NUM_SUB_EDGES * NUM_SUB_EDGES; }
+int MTriangle6::getNumFacesRep(){ return SQU(CTX.mesh.num_sub_edges); }
+int MTriangleN::getNumFacesRep(){ return SQU(CTX.mesh.num_sub_edges); }
-static void _myGetFaceRep(MTriangle *t, int num, double *x, double *y, double *z, SVector3 *n,
- int numSubEdges)
+static void _myGetFaceRep(MTriangle *t, int num, double *x, double *y, double *z,
+ SVector3 *n, int numSubEdges)
{
- // on the first layer, we have (numSubEdges-1) * 2 + 1 triangles
- // on the second layer, we have (numSubEdges-2) * 2 + 1 triangles
- // on the ith layer, we have (numSubEdges-1-i) * 2 + 1 triangles
+ // on the first layer, we have (numSubEdges-1) * 2 + 1 triangles
+ // on the second layer, we have (numSubEdges-2) * 2 + 1 triangles
+ // on the ith layer, we have (numSubEdges-1-i) * 2 + 1 triangles
int ix = 0, iy = 0;
int nbt = 0;
for (int i = 0; i < numSubEdges; i++){
@@ -926,14 +927,15 @@ static void _myGetFaceRep(MTriangle *t, int num, double *x, double *y, double *z
z[0] = pnt1.z(); z[1] = pnt2.z(); z[2] = pnt3.z();
}
-void MTriangleN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){
- _myGetFaceRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+void MTriangleN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+{
+ _myGetFaceRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
-void MTriangle6::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){
- _myGetFaceRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+void MTriangle6::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
+{
+ _myGetFaceRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
-
void MTriangle::getIntegrationPoints(int pOrder, int *npts, IntPt **pts) const
{
#if !defined(HAVE_GMSH_EMBEDDED)
@@ -980,7 +982,7 @@ double MTetrahedronN::distoShapeMeasure()
#if defined(HAVE_GMSH_EMBEDDED)
return 1.;
#else
- // if (_disto < -1.e21)
+ // if (_disto < -1.e21)
_disto = qmDistorsionOfMapping(this);
return _disto;
#endif
@@ -1052,10 +1054,11 @@ const gmshFunctionSpace* MTetrahedron::getFunctionSpace(int o) const
return 0;
}
-int MTetrahedron10::getNumEdgesRep(){ return 6 * NUM_SUB_EDGES; }
-int MTetrahedronN::getNumEdgesRep(){ return 6 * NUM_SUB_EDGES; }
+int MTetrahedron10::getNumEdgesRep(){ return 6 * CTX.mesh.num_sub_edges; }
+int MTetrahedronN::getNumEdgesRep(){ return 6 * CTX.mesh.num_sub_edges; }
-static void _myGetEdgeRep(MTetrahedron *tet, int num, double *x, double *y, double *z, SVector3 *n, int numSubEdges)
+static void _myGetEdgeRep(MTetrahedron *tet, int num, double *x, double *y, double *z,
+ SVector3 *n, int numSubEdges)
{
static double pp[4][3] = {{0,0,0},{1,0,0},{0,1,0},{0,0,1}};
static int ed [6][2] = {{0,1},{0,2},{0,3},{1,2},{1,3},{2,3}};
@@ -1080,31 +1083,36 @@ static void _myGetEdgeRep(MTetrahedron *tet, int num, double *x, double *y, doub
x[0] = pnt1.x(); x[1] = pnt2.x();
y[0] = pnt1.y(); y[1] = pnt2.y();
z[0] = pnt1.z(); z[1] = pnt2.z();
+
+ // not great, but better than nothing
+ static const int f[6] = {0, 0, 0, 1, 2, 3};
+ n[0] = n[1] = tet->getFace(f[iEdge]).normal();
}
void MTetrahedron10::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
{
- _myGetEdgeRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+ _myGetEdgeRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
void MTetrahedronN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
{
- _myGetEdgeRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+ _myGetEdgeRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
-int MTetrahedronN::getNumFacesRep(){ return 4 * NUM_SUB_EDGES * NUM_SUB_EDGES; }
-int MTetrahedron10::getNumFacesRep(){ return 4 * NUM_SUB_EDGES * NUM_SUB_EDGES; }
+int MTetrahedronN::getNumFacesRep(){ return 4 * SQU(CTX.mesh.num_sub_edges); }
+int MTetrahedron10::getNumFacesRep(){ return 4 * SQU(CTX.mesh.num_sub_edges); }
-void _myGetFaceRep (MTetrahedron *tet, int num, double *x, double *y, double *z, SVector3 *n, int numSubEdges )
+static void _myGetFaceRep(MTetrahedron *tet, int num, double *x, double *y, double *z,
+ SVector3 *n, int numSubEdges)
{
static double pp[4][3] = {{0,0,0},{1,0,0},{0,1,0},{0,0,1}};
static int fak [4][3] = {{0,1,2},{0,1,3},{0,2,3},{1,2,3}};
int iFace = num / (numSubEdges * numSubEdges);
int iSubFace = num % (numSubEdges * numSubEdges);
- int iVertex1 = fak [iFace][0];
- int iVertex2 = fak [iFace][1];
- int iVertex3 = fak [iFace][2];
+ int iVertex1 = fak[iFace][0];
+ int iVertex2 = fak[iFace][1];
+ int iVertex3 = fak[iFace][2];
/*
0
@@ -1115,9 +1123,9 @@ void _myGetFaceRep (MTetrahedron *tet, int num, double *x, double *y, double *z,
0 1 2 3 4 5
*/
- // on the first layer, we have (numSubEdges-1) * 2 + 1 triangles
- // on the second layer, we have (numSubEdges-2) * 2 + 1 triangles
- // on the ith layer, we have (numSubEdges-1-i) * 2 + 1 triangles
+ // on the first layer, we have (numSubEdges-1) * 2 + 1 triangles
+ // on the second layer, we have (numSubEdges-2) * 2 + 1 triangles
+ // on the ith layer, we have (numSubEdges-1-i) * 2 + 1 triangles
int ix = 0, iy = 0;
int nbt = 0;
for (int i = 0; i < numSubEdges; i++){
@@ -1166,45 +1174,22 @@ void _myGetFaceRep (MTetrahedron *tet, int num, double *x, double *y, double *z,
y[0] = pnt1.y(); y[1] = pnt2.y(); y[2] = pnt3.y();
z[0] = pnt1.z(); z[1] = pnt2.z(); z[2] = pnt3.z();
- // facetted first
-
SVector3 d1(x[1]-x[0],y[1]-y[0],z[1]-z[0]);
SVector3 d2(x[2]-x[0],y[2]-y[0],z[2]-z[0]);
n[0] = crossprod(d1, d2);
n[0].normalize();
n[1] = n[0];
n[2] = n[0];
-
- return;
-
- {
- SVector3 d1(J1[0][0], J1[0][1], J1[0][2]);
- SVector3 d2(J1[1][0], J1[1][1], J1[1][2]);
- n[0] = crossprod(d1, d2);
- n[0].normalize();
- }
- {
- SVector3 d1(J2[0][0], J2[0][1], J2[0][2]);
- SVector3 d2(J2[1][0], J2[1][1], J2[1][2]);
- n[1] = crossprod(d1, d2);
- n[1].normalize();
- }
- {
- SVector3 d1(J3[0][0], J3[0][1], J3[0][2]);
- SVector3 d2(J3[1][0], J3[1][1], J3[1][2]);
- n[2] = crossprod(d1, d2);
- n[2].normalize();
- }
}
void MTetrahedronN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
{
- _myGetFaceRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+ _myGetFaceRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
void MTetrahedron10::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
{
- _myGetFaceRep(this,num,x,y,z,n,NUM_SUB_EDGES);
+ _myGetFaceRep(this, num, x, y, z, n, CTX.mesh.num_sub_edges);
}
void MTetrahedron::getIntegrationPoints(int pOrder, int *npts, IntPt **pts) const
diff --git a/doc/texinfo/gmsh.texi b/doc/texinfo/gmsh.texi
index 7880e86195..46479213dc 100644
--- a/doc/texinfo/gmsh.texi
+++ b/doc/texinfo/gmsh.texi
@@ -1999,12 +1999,18 @@ characteristic lengths:
@ftable @code
@item Characteristic Length @{ @var{expression-list} @} = @var{expression};
-Modifies the characteristic length of the points whose identification
+Modify the characteristic length of the points whose identification
numbers are listed in @var{expression-list}. The new value is given by
@var{expression}.
@item Field[@var{expression}] = @var{string};
+Create a new field (with id number @var{expression}), of type
+@var{string}.
@item Field[@var{expression}].@var{string} = @var{char-expression} | @var{expression} | @var{expression-list};
+Set the option @var{string} of the @var{expression}-th field.
@item Background Field = @var{expression};
+Select the @var{expression}-th field as the one used to compute element
+sizes. Only one background field can be given; if you want to combine
+several field, use the @code{Min} or @code{Max} field (see below).
@end ftable
Here is the list of all available fields with their associated options:
@@ -2088,7 +2094,7 @@ line.
(A deprecated synonym for @code{Progression} is @code{Power}.)
-@item Transfinite Surface @{ @var{expression} @} = @{ @var{expression-list} @} < Left | Right | Alternate > ;
+@item Transfinite Surface @{ @var{expression} @} < = @{ @var{expression-list} @} > < Left | Right | Alternate > ;
Selects the surface @var{expression} to be meshed with the 2D
transfinite algorithm. The @var{expression-list} should contain the
identification numbers of three or four points on the boundary of the
--
GitLab