diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index fb67f3f7cf31aadbc448a97f8f5f25cd952e748f..d5f0c88040b28a181767f4813204bc5439b3d364 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -8,12 +8,12 @@
#include "MElement.h"
#include "GEntity.h"
#include "GFace.h"
+#include "FunctionSpace.h"
#if defined(HAVE_GMSH_EMBEDDED)
# include "GmshEmbedded.h"
#else
# include "Numeric.h"
-# include "FunctionSpace.h"
# include "GaussLegendre1D.h"
# include "Message.h"
# include "Context.h"
@@ -423,25 +423,20 @@ double MElement::interpolateDiv(double val[], double u, double v, double w, int
interpolateGrad(&val[2], u, v, w, fz, stride, inv);
return fx[0] + fy[1] + fz[2];
}
-
-
-// Expensive, generic version
-
-void MElement::pnt(double u,double v,double w,SPoint3& p) {
- double x(0),y(0),z(0);
-
- for (int i=0;i<getNumVertices();i++) {
+// Expensive, generic version
+void MElement::pnt(double u, double v, double w, SPoint3 &p)
+{
+ double x(0), y(0), z(0);
+ for (int i = 0; i < getNumVertices(); i++) {
double s;
- getShapeFunction(i,u,v,w,s);
-
+ getShapeFunction(i, u, v, w, s);
MVertex* v = getVertex(i);
-
x += v->x() * s;
y += v->y() * s;
z += v->z() * s;
}
- p = SPoint3(x,y,z);
+ p = SPoint3(x, y, z);
}
void MElement::writeMSH(FILE *fp, double version, bool binary, int num,
@@ -742,20 +737,21 @@ int MElement::getInfoMSH(const int typeMSH, const char **const name)
}
}
-void MLine::pnt(double u,double v,double w,SPoint3& p) {
- double f1,f2;
- getShapeFunction(0,u,v,w,f1);
- getShapeFunction(1,u,v,w,f2);
+void MLine::pnt(double u,double v,double w,SPoint3& p)
+{
+ double f1, f2;
+ getShapeFunction(0, u, v, w, f1);
+ getShapeFunction(1, u, v, w, f2);
double x = f1 * _v[0]->x() + f2 * _v[1]->x();
double y = f1 * _v[0]->y() + f2 * _v[1]->y();
double z = f1 * _v[0]->z() + f2 * _v[1]->z();
- p = SPoint3(x,y,z);
+ p = SPoint3(x, y, z);
}
-void MLine3::pnt(double u, double v, double w, SPoint3 &p) {
-
+void MLine3::pnt(double u, double v, double w, SPoint3 &p)
+{
double sf[3];
gmshFunctionSpaces::find(MSH_LIN_3).f(u, v, w, sf);
double x = sf[0] * _v[0]->x() + sf[1] * _v[1]->x() + sf[2] * _vs[0]->x();
@@ -765,35 +761,24 @@ void MLine3::pnt(double u, double v, double w, SPoint3 &p) {
p = SPoint3(x,y,z);
}
-void MLine3::getShapeFunction(int num,double uu, double vv,double ww,double& s) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
-
+void MLine3::getShapeFunction(int num,double uu, double vv,double ww,double& s)
+{
if (num > 2) s = 0.;
double sf[3];
gmshFunctionSpaces::find(MSH_LIN_3).f(uu, vv, ww, sf);
s = sf[num];
-
-#endif
-
}
-void MLine3::getGradShapeFunction(int num,double uu,double vv,double ww,double s[3]) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
-
- if (num > 2) for (int i=0;i<3;i++) s[i] = 0.;
+void MLine3::getGradShapeFunction(int num,double uu,double vv,double ww,double s[3])
+{
+ if (num > 2) for (int i = 0; i < 3; i++) s[i] = 0.;
double sf[3][3];
gmshFunctionSpaces::find(MSH_LIN_3).df(uu, vv, ww, sf);
- for (int i=0;i<3;i++) s[i] = sf[num][i];
-
-#endif
-
+ for (int i = 0; i < 3; i++) s[i] = sf[num][i];
}
-
-void MLineN::pnt(double uu, double vv, double ww, SPoint3 &p) {
-
+void MLineN::pnt(double uu, double vv, double ww, SPoint3 &p)
+{
double sf[6];
switch (getPolynomialOrder()) {
@@ -802,25 +787,25 @@ void MLineN::pnt(double uu, double vv, double ww, SPoint3 &p) {
case 3: gmshFunctionSpaces::find(MSH_LIN_4).f(uu, vv, ww, sf); break;
case 4: gmshFunctionSpaces::find(MSH_LIN_5).f(uu, vv, ww, sf); break;
case 5: gmshFunctionSpaces::find(MSH_LIN_6).f(uu, vv, ww, sf); break;
- default: Msg::Error("Order %d line point interpolation not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d line point interpolation not implemented",
+ getPolynomialOrder());
+ break;
}
- double x = 0; for (int i=0;i<2;i++) x += sf[i] * _v[i]->x();
- double y = 0; for (int i=0;i<2;i++) y += sf[i] * _v[i]->y();
- double z = 0; for (int i=0;i<2;i++) z += sf[i] * _v[i]->z();
+ double x = 0; for (int i = 0; i < 2; i++) x += sf[i] * _v[i]->x();
+ double y = 0; for (int i = 0; i < 2; i++) y += sf[i] * _v[i]->y();
+ double z = 0; for (int i = 0; i < 2; i++) z += sf[i] * _v[i]->z();
- for (int i=0;i<_vs.size();i++) x += sf[i+2] * _vs[i]->x();
- for (int i=0;i<_vs.size();i++) y += sf[i+2] * _vs[i]->y();
- for (int i=0;i<_vs.size();i++) z += sf[i+2] * _vs[i]->z();
+ for (int i = 0; i < _vs.size(); i++) x += sf[i+2] * _vs[i]->x();
+ for (int i = 0; i < _vs.size(); i++) y += sf[i+2] * _vs[i]->y();
+ for (int i = 0; i < _vs.size(); i++) z += sf[i+2] * _vs[i]->z();
- p = SPoint3(x,y,z);
-
+ p = SPoint3(x, y, z);
}
-void MLineN::getShapeFunction(int num,double uu, double vv,double ww,double& s) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
-
+void MLineN::getShapeFunction(int num, double uu, double vv, double ww, double &s)
+{
double sf[6];
switch (getPolynomialOrder()) {
case 1: gmshFunctionSpaces::find(MSH_LIN_2).f(uu, vv, ww, sf); break;
@@ -828,18 +813,16 @@ void MLineN::getShapeFunction(int num,double uu, double vv,double ww,double& s)
case 3: gmshFunctionSpaces::find(MSH_LIN_4).f(uu, vv, ww, sf); break;
case 4: gmshFunctionSpaces::find(MSH_LIN_5).f(uu, vv, ww, sf); break;
case 5: gmshFunctionSpaces::find(MSH_LIN_6).f(uu, vv, ww, sf); break;
- default: Msg::Error("Order %d line shape function not provided", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d line shape function not provided",
+ getPolynomialOrder());
+ break;
}
s = sf[num];
-
-#endif
-
}
-void MLineN::getGradShapeFunction(int num,double uu,double vv,double ww,double s[3]) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
-
+void MLineN::getGradShapeFunction(int num, double uu, double vv, double ww, double s[3])
+{
double sf[6][3];
switch (getPolynomialOrder()) {
case 1: gmshFunctionSpaces::find(MSH_LIN_2).df(uu, vv, ww, sf); break;
@@ -847,17 +830,16 @@ void MLineN::getGradShapeFunction(int num,double uu,double vv,double ww,double s
case 3: gmshFunctionSpaces::find(MSH_LIN_4).df(uu, vv, ww, sf); break;
case 4: gmshFunctionSpaces::find(MSH_LIN_5).df(uu, vv, ww, sf); break;
case 5: gmshFunctionSpaces::find(MSH_LIN_6).df(uu, vv, ww, sf); break;
- default: Msg::Error("Order %d line shape function gradient not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d line shape function gradient not implemented",
+ getPolynomialOrder());
+ break;
}
- for (int i=0;i<3;i++) s[i] = sf[num][i];
-#endif
-
+ for (int i = 0; i < 3; i++) s[i] = sf[num][i];
}
-void MTriangle::getShapeFunction(int num,double uu, double vv,double ww,double& s) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
-
+void MTriangle::getShapeFunction(int num, double uu, double vv, double ww, double &s)
+{
int nf = getNumFaceVertices();
double fv[256];
@@ -868,7 +850,10 @@ void MTriangle::getShapeFunction(int num,double uu, double vv,double ww,double&
case 3: gmshFunctionSpaces::find(MSH_TRI_9).f(uu, vv, fv); break;
case 4: gmshFunctionSpaces::find(MSH_TRI_12).f(uu, vv, fv); break;
case 5: gmshFunctionSpaces::find(MSH_TRI_15I).f(uu, vv, fv); break;
- default: Msg::Error("Order %d triangle shape function not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d triangle shape function not implemented",
+ getPolynomialOrder());
+ break;
}
}
else{
@@ -878,18 +863,17 @@ void MTriangle::getShapeFunction(int num,double uu, double vv,double ww,double&
case 3: gmshFunctionSpaces::find(MSH_TRI_10).f(uu, vv, fv); break;
case 4: gmshFunctionSpaces::find(MSH_TRI_15).f(uu, vv, fv); break;
case 5: gmshFunctionSpaces::find(MSH_TRI_21).f(uu, vv, fv); break;
- default: Msg::Error("Order %d triangle shape function not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d triangle shape function not implemented",
+ getPolynomialOrder());
+ break;
}
}
s = fv[num];
-#endif
-
}
-void MTriangle::getGradShapeFunction(int num,double uu,double vv,double ww,double s[3]) {
-
-
-#if !defined(HAVE_GMSH_EMBEDDED)
+void MTriangle::getGradShapeFunction(int num, double uu, double vv, double ww, double s[3])
+{
double grads[256][2];
int nf = getNumFaceVertices();
@@ -900,7 +884,10 @@ void MTriangle::getGradShapeFunction(int num,double uu,double vv,double ww,doubl
case 3: gmshFunctionSpaces::find(MSH_TRI_9).df(uu, vv, ww, grads); break;
case 4: gmshFunctionSpaces::find(MSH_TRI_12).df(uu, vv, ww, grads); break;
case 5: gmshFunctionSpaces::find(MSH_TRI_15I).df(uu, vv, ww, grads); break;
- default: Msg::Error("Order %d triangle shape function gradient not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d triangle shape function gradient not implemented",
+ getPolynomialOrder());
+ break;
}
}
else{
@@ -910,19 +897,20 @@ void MTriangle::getGradShapeFunction(int num,double uu,double vv,double ww,doubl
case 3: gmshFunctionSpaces::find(MSH_TRI_10).df(uu, vv, ww, grads); break;
case 4: gmshFunctionSpaces::find(MSH_TRI_15).df(uu, vv, ww, grads); break;
case 5: gmshFunctionSpaces::find(MSH_TRI_21).df(uu, vv, ww, grads); break;
- default: Msg::Error("Order %d triangle shape function gradient not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d triangle shape function gradient not implemented",
+ getPolynomialOrder());
+ break;
}
}
- for (int i=0;i<2;i++) s[i] = grads[num][i];
+ for (int i = 0; i < 2; i++) s[i] = grads[num][i];
s[2] = 0;
-#endif
}
void MTriangle::jac(int ord, MVertex *vs[], double uu, double vv, double ww,
double j[2][3])
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double grads[256][2];
int nf = getNumFaceVertices();
@@ -961,13 +949,11 @@ void MTriangle::jac(int ord, MVertex *vs[], double uu, double vv, double ww,
for(int i = 3; i < 3 * ord + nf; i++) j[1][1] += grads[i][1] * vs[i - 3]->y();
for(int i = 3; i < 3 * ord + nf; i++) j[0][2] += grads[i][0] * vs[i - 3]->z();
for(int i = 3; i < 3 * ord + nf; i++) j[1][2] += grads[i][1] * vs[i - 3]->z();
-#endif
}
void MTriangle::pnt(int ord, MVertex *vs[], double uu, double vv, double ww,
SPoint3 &p)
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double sf[256];
int nf = getNumFaceVertices();
@@ -1000,14 +986,12 @@ void MTriangle::pnt(int ord, MVertex *vs[], double uu, double vv, double ww,
for(int i = 3; i < 3 * ord + nf; i++) y += sf[i] * vs[i - 3]->y();
for(int i = 3; i < 3 * ord + nf; i++) z += sf[i] * vs[i - 3]->z();
- p = SPoint3(x,y,z);
-#endif
+ p = SPoint3(x, y, z);
}
void MTetrahedron::jac(int ord, MVertex *vs[], double uu, double vv, double ww,
double j[3][3])
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double grads[256][3];
switch(ord){
case 1: gmshFunctionSpaces::find(MSH_TET_4) .df(uu, vv, ww, grads); break;
@@ -1045,13 +1029,11 @@ void MTetrahedron::jac(int ord, MVertex *vs[], double uu, double vv, double ww,
for(int i = 4; i < N; i++) j[0][2] += grads[i][0] * vs[i - 4]->z();
for(int i = 4; i < N; i++) j[1][2] += grads[i][1] * vs[i - 4]->z();
for(int i = 4; i < N; i++) j[2][2] += grads[i][2] * vs[i - 4]->z();
-#endif
}
void MTetrahedron::pnt(int ord, MVertex *vs[], double uu, double vv, double ww,
SPoint3 &p)
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double sf[256];
int nv = getNumVolumeVertices();
@@ -1078,8 +1060,8 @@ void MTetrahedron::pnt(int ord, MVertex *vs[], double uu, double vv, double ww,
double y = 0 ; for(int i = 0; i < 4; i++) y += sf[i] * _v[i]->y();
double z = 0 ; for(int i = 0; i < 4; i++) z += sf[i] * _v[i]->z();
- int nbComplete = (ord+1)*(ord+2)*(ord+3)/6;
- int nbInterior = (ord-3)*(ord-2)*(ord-1)/6;
+ int nbComplete = (ord + 1) * (ord + 2) * (ord + 3) / 6;
+ int nbInterior = (ord - 3) * (ord - 2) * (ord - 1) / 6;
const int N = nv ? nbComplete : nbComplete - nbInterior;
@@ -1087,8 +1069,7 @@ void MTetrahedron::pnt(int ord, MVertex *vs[], double uu, double vv, double ww,
for(int i = 4; i < N; i++) y += sf[i] * vs[i - 4]->y();
for(int i = 4; i < N; i++) z += sf[i] * vs[i - 4]->z();
- p = SPoint3(x,y,z);
-#endif
+ p = SPoint3(x, y, z);
}
const int numSubEdges = 6;
@@ -1192,28 +1173,23 @@ void MTriangleN::getEdgeRep(int num, double *x, double *y, double *z, SVector3 *
}
}
-void MTetrahedron10::getShapeFunction(int num, double uu, double vv, double ww, double &s) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
+void MTetrahedron10::getShapeFunction(int num, double uu, double vv, double ww, double &s)
+{
double sf[256];
gmshFunctionSpaces::find(MSH_TET_10).f(uu, vv, ww, sf);
s = sf[num];
-#endif
}
-void MTetrahedron10::getGradShapeFunction(int num, double uu, double vv, double ww, double gs[3]) {
-
-#if !defined(HAVE_GMSH_EMBEDDED)
+void MTetrahedron10::getGradShapeFunction(int num, double uu, double vv, double ww,
+ double gs[3])
+{
double gsf[256][3];
gmshFunctionSpaces::find(MSH_TET_10).df(uu, vv, ww, gsf);
- for (int i=0;i<3;i++) gs[i] = gsf[num][i];
-#endif
+ for (int i = 0; i < 3; i++) gs[i] = gsf[num][i];
}
-
void MTetrahedronN::getShapeFunction(int num, double uu, double vv, double ww, double &s)
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double sf[256];
int nv = getNumVolumeVertices();
@@ -1225,24 +1201,27 @@ void MTetrahedronN::getShapeFunction(int num, double uu, double vv, double ww, d
case 3: gmshFunctionSpaces::find(MSH_TET_20).f(uu, vv, ww, sf); break;
case 4: gmshFunctionSpaces::find(MSH_TET_34).f(uu, vv, ww, sf); break;
case 5: gmshFunctionSpaces::find(MSH_TET_52).f(uu, vv, ww, sf); break;
- default: Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder());
+ break;
}
}
else{
switch(getPolynomialOrder()){
case 4: gmshFunctionSpaces::find(MSH_TET_35).f(uu, vv, ww, sf); break;
case 5: gmshFunctionSpaces::find(MSH_TET_56).f(uu, vv, ww, sf); break;
- default: Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder());
+ break;
}
}
s = sf[num];
-#endif
}
-void MTetrahedronN::getGradShapeFunction(int num, double uu, double vv, double ww, double gs[3])
+void MTetrahedronN::getGradShapeFunction(int num, double uu, double vv, double ww,
+ double gs[3])
{
-#if !defined(HAVE_GMSH_EMBEDDED)
double gsf[256][3];
int nv = getNumVolumeVertices();
@@ -1254,18 +1233,21 @@ void MTetrahedronN::getGradShapeFunction(int num, double uu, double vv, double w
case 3: gmshFunctionSpaces::find(MSH_TET_20).df(uu, vv, ww, gsf); break;
case 4: gmshFunctionSpaces::find(MSH_TET_34).df(uu, vv, ww, gsf); break;
case 5: gmshFunctionSpaces::find(MSH_TET_52).df(uu, vv, ww, gsf); break;
- default: Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder());
+ break;
}
}
else{
switch(getPolynomialOrder()){
case 4: gmshFunctionSpaces::find(MSH_TET_35).df(uu, vv, ww, gsf); break;
case 5: gmshFunctionSpaces::find(MSH_TET_56).df(uu, vv, ww, gsf); break;
- default: Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder()); break;
+ default:
+ Msg::Error("Order %d tetrahedron pnt not implemented", getPolynomialOrder());
+ break;
}
}
- for (int i=0;i<3;i++) gs[i] = gsf[num][i];
-#endif
+ for (int i = 0; i < 3; i++) gs[i] = gsf[num][i];
}
int MTetrahedronN::getNumEdgesRep(){ return 6 * numSubEdges; }
@@ -1298,7 +1280,6 @@ void MTetrahedronN::getEdgeRep(int num, double *x, double *y, double *z, SVector
z[0] = pnt1.z(); z[1] = pnt2.z();
}
-
int MTetrahedronN::getNumFacesRep(){ return 4 * numSubEdges * numSubEdges ; }
void MTetrahedronN::getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
diff --git a/Geo/MElement.h b/Geo/MElement.h
index 060c4189f30edb7850e52f92417f95775e676abe..1c5d2111c1d0bb1ec3960c92a6f36782835025b2 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -15,7 +15,6 @@
#include "MFace.h"
#include "Message.h"
-
struct IntPt{
double pt[3];
double weight;
@@ -342,10 +341,6 @@ class MLine : public MElement {
{
MVertex *tmp = _v[0]; _v[0] = _v[1]; _v[1] = tmp;
}
-
-
- virtual void pnt(double u, double v, double w, SPoint3 &p);
-
virtual void getShapeFunction(int num, double u, double v, double w, double &s)
{
switch(num) {
@@ -362,6 +357,7 @@ class MLine : public MElement {
default : s[0] = s[1] = s[2] = 0.; break;
}
}
+ virtual void pnt(double u, double v, double w, SPoint3 &p);
virtual bool isInside(double u, double v, double w, double tol=1.e-8)
{
if(u < -(1. + tol) || u > (1. + tol))
@@ -418,12 +414,9 @@ class MLine3 : public MLine {
MLine::_getEdgeVertices(v);
v[2] = _vs[0];
}
-
- virtual void pnt(double u, double v, double w, SPoint3 &p);
-
virtual void getShapeFunction(int num, double u, double v, double w, double &s);
virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]);
-
+ virtual void pnt(double u, double v, double w, SPoint3 &p);
virtual int getTypeForMSH(){ return MSH_LIN_3; }
virtual int getTypeForUNV(){ return 24; } // parabolic beam
virtual int getTypeForVTK(){ return 21; }
@@ -472,12 +465,9 @@ class MLineN : public MLine {
MLine::_getEdgeVertices(v);
for(unsigned int i = 0; i != _vs.size(); ++i) v[i+2] = _vs[i];
}
-
- virtual void pnt(double u, double v, double w, SPoint3 &p) ;
-
virtual void getShapeFunction(int num, double u, double v, double w, double &s);
virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]);
-
+ virtual void pnt(double u, double v, double w, SPoint3 &p);
virtual int getTypeForMSH()
{
if(_vs.size() == 2) return MSH_LIN_4;
@@ -585,34 +575,8 @@ class MTriangle : public MElement {
{
MVertex *tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
}
-
virtual void getShapeFunction(int num, double u, double v, double w, double &s);
virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]);
-
-/* virtual void getShapeFunction(int num, double u, double v, double w, double &s) */
-/* { */
-/* switch(num){ */
-/* case 0 : s = 1. - u - v; break; */
-/* case 1 : s = u ; break; */
-/* case 2 : s = v; break; */
-/* default : s = 0.; break; */
-/* } */
-/* } */
-/* virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]) */
-/* { */
-/* switch(num) { */
-/* case 0 : s[0] = -1.; s[1] = -1.; s[2] = 0.; break; */
-/* case 1 : s[0] = 1.; s[1] = 0.; s[2] = 0.; break; */
-/* case 2 : s[0] = 0.; s[1] = 1.; s[2] = 0.; break; */
-/* default : s[0] = s[1] = s[2] = 0.; break; */
-/* } */
-/* } */
- virtual bool isInside(double u, double v, double w, double tol=1.e-8)
- {
- if(u < (-tol) || v < (-tol) || u > ((1. + tol) - v))
- return false;
- return true;
- }
virtual void jac(int order, MVertex *verts[], double u, double v, double w,
double j[2][3]);
virtual void jac(double u, double v, double w, double j[2][3])
@@ -625,6 +589,12 @@ class MTriangle : public MElement {
{
pnt(1, 0, u, v, w, p);
}
+ virtual bool isInside(double u, double v, double w, double tol=1.e-8)
+ {
+ if(u < (-tol) || v < (-tol) || u > ((1. + tol) - v))
+ return false;
+ return true;
+ }
virtual void getIntegrationPoints(int pOrder, int *npts, IntPt **pts) const;
virtual SPoint3 circumcenter();
private:
@@ -727,7 +697,6 @@ class MTriangle6 : public MTriangle {
tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
tmp = _vs[0]; _vs[0] = _vs[2]; _vs[2] = tmp;
}
-
virtual void jac(double u, double v, double w, double j[2][3])
{
MTriangle::jac(2, _vs, u, v, w, j);
@@ -1506,9 +1475,6 @@ class MTetrahedron10 : public MTetrahedron {
*
*/
-
-
-
static int reverseTet20[20] = {0,2,1,3, // principal vertices
9,8, // E0 switches with E2
7,6, // E1 inverts direction
@@ -1552,14 +1518,14 @@ class MTetrahedronN : public MTetrahedron {
std::vector<MVertex *> _vs;
const short _order;
public:
- MTetrahedronN(MVertex *v0, MVertex *v1, MVertex *v2, MVertex *v3,
- std::vector<MVertex*> &v, int order, int num=0, int part=0)
- : MTetrahedron(v0, v1, v2, v3, num, part) , _vs (v), _order(order)
+ MTetrahedronN(MVertex *v0, MVertex *v1, MVertex *v2, MVertex *v3,
+ std::vector<MVertex*> &v, int order, int num=0, int part=0)
+ : MTetrahedron(v0, v1, v2, v3, num, part) , _vs (v), _order(order)
{
for(unsigned int i = 0; i < _vs.size(); i++) _vs[i]->setPolynomialOrder(_order);
}
- MTetrahedronN(std::vector<MVertex*> &v, int order, int num=0, int part=0)
- : MTetrahedron(v[0], v[1], v[2], v[3], num, part) , _order(order)
+ MTetrahedronN(std::vector<MVertex*> &v, int order, int num=0, int part=0)
+ : MTetrahedron(v[0], v[1], v[2], v[3], num, part) , _order(order)
{
for(unsigned int i = 4; i < v.size(); i++) _vs.push_back(v[i]);
for(unsigned int i = 0; i < _vs.size(); i++) _vs[i]->setPolynomialOrder(_order);
@@ -1571,15 +1537,15 @@ class MTetrahedronN : public MTetrahedron {
virtual int getNumEdgeVertices() const { return _order - 1; }
virtual int getNumFaceVertices() const
{
- return ((_order - 1)*(_order - 2))/2;
+ return ((_order - 1) * (_order - 2)) / 2;
}
virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(2 + getNumEdgeVertices());
MTetrahedron::_getEdgeVertices(num, v);
int j = 2;
- const int ie = (num+1)*getNumEdgeVertices();
- for(int i = num*getNumEdgeVertices(); i != ie; ++i) v[j++] = _vs[i];
+ const int ie = (num + 1) * getNumEdgeVertices();
+ for(int i = num * getNumEdgeVertices(); i != ie; ++i) v[j++] = _vs[i];
}
virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
{
@@ -1598,7 +1564,6 @@ class MTetrahedronN : public MTetrahedron {
}
}
virtual int getNumEdgeVertices(){ return _order - 1; }
-
virtual int getTypeForMSH()
{
// (p+1)*(p+2)*(p+3)/6
@@ -1643,15 +1608,12 @@ class MTetrahedronN : public MTetrahedron {
}
}
}
-
virtual void getShapeFunction(int num, double u, double v, double w, double &s);
virtual void getGradShapeFunction(int num, double u, double v, double w, double s[3]);
-
virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual int getNumEdgesRep();
virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual int getNumFacesRep();
-
virtual void jac(double u, double v, double w, double j[3][3])
{
MTetrahedron::jac(_order, _vs , u, v, w, j);
diff --git a/Makefile b/Makefile
index e7f9bcab69b41b92297e475ea291d6bedf347882..88d66dc695c1e0207002100da0d49ac5f4486dda 100644
--- a/Makefile
+++ b/Makefile
@@ -32,7 +32,8 @@ GMSH_EMBEDDED = ${GMSH_API} Geo/discrete*.cpp\
Geo/GEntity.cpp Geo/GVertex.cpp Geo/GEdge.cpp\
Geo/GEdgeLoop.cpp Geo/GFace.cpp Geo/GRegion.cpp\
Geo/MElement.cpp Geo/MFace.cpp Geo/MVertex.cpp\
- Common/StringUtils.{cpp,h} Numeric/NumericEmbedded.{cpp,h}\
+ Common/StringUtils.{cpp,h}\
+ Numeric/NumericEmbedded.{cpp,h} Numeric/FunctionSpace.{cpp,h}\
utils/embed/GmshEmbedded.{cpp,h} utils/embed/Makefile
all: link
diff --git a/Mesh/BackgroundMesh.cpp b/Mesh/BackgroundMesh.cpp
index ee68065662d849b8f05b7f43c557d0a822d84ded..959f548ff21ae51148edab09064e75ac9c418f31 100644
--- a/Mesh/BackgroundMesh.cpp
+++ b/Mesh/BackgroundMesh.cpp
@@ -162,10 +162,9 @@ double BGM_MeshSize(GEntity *ge, double U, double V, double X, double Y, double
lc = std::min(lc, CTX.mesh.lc_max);
if(lc <= 0.){
- Msg::Error("Wrong characteristic length lc = %g", lc);
- printf("%g %g \n", CTX.mesh.lc_min, CTX.mesh.lc_max);
-
- lc = l1;
+ Msg::Error("Wrong characteristic length lc = %g (lcmin = %g, lcmax = %g)",
+ lc, CTX.mesh.lc_min, CTX.mesh.lc_max);
+ lc = l1;
}
return lc * CTX.mesh.lc_factor;
diff --git a/Mesh/Field.h b/Mesh/Field.h
index c6333532ad11896ef46bd12d7650314d496069ed..00abcd052edf5e86765350bad4618fac5ab30ca0 100644
--- a/Mesh/Field.h
+++ b/Mesh/Field.h
@@ -49,7 +49,7 @@ class FieldOption {
case FIELD_OPTION_STRING: typ = "string"; break;
case FIELD_OPTION_LIST: typ = "list"; break;
}
- return std::string(_help + " (type: " + typ + "; default value: " + val + ")");
+ return _help + " (type: " + typ + "; default value: " + val + ")";
}
virtual void numerical_value(double val) { throw(1); }
virtual double numerical_value() const { throw(1); }
diff --git a/Numeric/FunctionSpace.cpp b/Numeric/FunctionSpace.cpp
index 76640414c32681fb37949dbdbd045fd3d324263d..85647e93e8dd1d59ffddfed14385d9acbedb477f 100644
--- a/Numeric/FunctionSpace.cpp
+++ b/Numeric/FunctionSpace.cpp
@@ -7,24 +7,22 @@
#include "GmshDefines.h"
#include "Message.h"
-Double_Matrix generate1DMonomials(int order) {
-
- Double_Matrix monomials(order+1,1);
- for (int i=0;i<order+1;i++) monomials(i,0) = i;
-
+Double_Matrix generate1DMonomials(int order)
+{
+ Double_Matrix monomials(order + 1, 1);
+ for (int i = 0; i < order + 1; i++) monomials(i, 0) = i;
+ return monomials;
}
-Double_Matrix generate1DPoints(int order) {
-
- Double_Matrix line(order+1,1);
-
- line(0,0) = -1.;
- line(1,0) = 1.;
-
- double dd = 2./order;
-
- for (int i=2;i<order;i++) line(i,0) = -1. + dd * (i-1);
-}
+Double_Matrix generate1DPoints(int order)
+{
+ Double_Matrix line(order + 1, 1);
+ line(0, 0) = -1.;
+ line(1, 0) = 1.;
+ double dd = 2. / order;
+ for (int i = 2; i < order; i++) line(i, 0) = -1. + dd * (i - 1);
+ return line;
+}
Double_Matrix generatePascalTriangle(int order)
{
@@ -166,7 +164,7 @@ int nbdoftriangleserendip(int order) { return 3 * order; }
// to numbering of principal vertices of face !!!!
// uv surface - orientation v0-v2-v1
-void nodepositionface0(int order, double *u, double *v, double *w)
+void nodepositionface0(int order, double *u, double *v, double *w)
{
int ndofT = nbdoftriangle(order);
if (order == 0) { u[0] = 0.; v[0] = 0.; w[0] = 0.; return; }
@@ -206,8 +204,9 @@ void nodepositionface0(int order, double *u, double *v, double *w)
w[k] = w[k] / order;
}
}
+
// uw surface - orientation v0-v1-v3
-void nodepositionface1(int order, double *u, double *v, double *w)
+void nodepositionface1(int order, double *u, double *v, double *w)
{
int ndofT = nbdoftriangle(order);
if (order == 0) { u[0] = 0.; v[0] = 0.; w[0] = 0.; return; }
@@ -247,7 +246,7 @@ void nodepositionface1(int order, double *u, double *v, double *w)
}
// vw surface - orientation v0-v3-v2
-void nodepositionface2(int order, double *u, double *v, double *w)
+void nodepositionface2(int order, double *u, double *v, double *w)
{
int ndofT = nbdoftriangle(order);
if (order == 0) { u[0] = 0.; v[0] = 0.; return; }
@@ -580,7 +579,7 @@ std::map<int, gmshFunctionSpace> gmshFunctionSpaces::fs;
const gmshFunctionSpace &gmshFunctionSpaces::find(int tag)
{
- std::map<int,gmshFunctionSpace>::const_iterator it = fs.find(tag);
+ std::map<int, gmshFunctionSpace>::const_iterator it = fs.find(tag);
if (it != fs.end()) return it->second;
gmshFunctionSpace F;