// Gmsh - Copyright (C) 1997-2010 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
#ifndef _MTRIANGLE_H_
#define _MTRIANGLE_H_
#include "MElement.h"
/*
* MTriangle
*
* v
* ^
* |
* 2
* |`\
* | `\
* | `\
* | `\
* | `\
* 0----------1 --> u
*
*/
class MTriangle : public MElement {
protected:
MVertex *_v[3];
void _getEdgeVertices(const int num, std::vector<MVertex*> &v) const
{
v[0] = _v[edges_tri(num, 0)];
v[1] = _v[edges_tri(num, 1)];
}
void _getFaceVertices(std::vector<MVertex*> &v) const
{
v[0] = _v[0];
v[1] = _v[1];
v[2] = _v[2];
}
public :
MTriangle(MVertex *v0, MVertex *v1, MVertex *v2, int num=0, int part=0)
: MElement(num, part)
{
_v[0] = v0; _v[1] = v1; _v[2] = v2;
}
MTriangle(std::vector<MVertex*> &v, int num=0, int part=0)
: MElement(num, part)
{
for(int i = 0; i < 3; i++) _v[i] = v[i];
}
~MTriangle(){}
virtual int getDim() const { return 2; }
virtual double etaShapeMeasure();
virtual double gammaShapeMeasure();
virtual double distoShapeMeasure();
virtual double getInnerRadius();
virtual double angleShapeMeasure();
virtual int getNumVertices() const { return 3; }
virtual MVertex *getVertex(int num){ return _v[num]; }
virtual MVertex *getVertexMED(int num)
{
static const int map[3] = {0, 2, 1};
return getVertex(map[num]);
}
virtual MVertex *getOtherVertex(MVertex *v1, MVertex *v2)
{
if(_v[0] != v1 && _v[0] != v2) return _v[0];
if(_v[1] != v1 && _v[1] != v2) return _v[1];
if(_v[2] != v1 && _v[2] != v2) return _v[2];
return 0;
}
virtual int getNumEdges(){ return 3; }
virtual MEdge getEdge(int num)
{
return MEdge(_v[edges_tri(num, 0)], _v[edges_tri(num, 1)]);
}
virtual void getEdgeInfo (const MEdge & edge, int &ithEdge, int &sign) const
{
for (ithEdge = 0; ithEdge < 3; ithEdge++){
const MVertex *v0 = _v[edges_tri(ithEdge, 0)];
const MVertex *v1 = _v[edges_tri(ithEdge, 1)];
if (v0 == edge.getVertex(0) && v1 == edge.getVertex(1)){
sign = 1; return;
}
if (v1 == edge.getVertex(0) && v0 == edge.getVertex(1)){
sign = -1; return;
}
}
Msg::Error("Could not get edge information for triangle %d", getNum());
}
virtual int getNumEdgesRep(){ return 3; }
virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
{
MEdge e(getEdge(num));
_getEdgeRep(e.getVertex(0), e.getVertex(1), x, y, z, n, 0);
}
virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(2);
_getEdgeVertices(num, v);
}
virtual int getNumFaces(){ return 1; }
virtual MFace getFace(int num)
{
return MFace(_v[0], _v[1], _v[2]);
}
virtual int getNumFacesRep(){ return 1; }
virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n)
{
_getFaceRep(_v[0], _v[1], _v[2], x, y, z, n);
}
virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(3);
_getFaceVertices(v);
}
virtual int getType() const { return TYPE_TRI; }
virtual int getTypeForMSH() const { return MSH_TRI_3; }
virtual int getTypeForUNV() const { return 91; } // thin shell linear triangle
virtual int getTypeForVTK() const { return 5; }
virtual const char *getStringForPOS() const { return "ST"; }
virtual const char *getStringForBDF() const { return "CTRIA3"; }
virtual const char *getStringForDIFF() const { return "ElmT3n2D"; }
virtual const char *getStringForINP() const { return "C2D3"; }
virtual void revert()
{
MVertex *tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
}
virtual const polynomialBasis* getFunctionSpace(int o=-1) const;
virtual const JacobianBasis* getJacobianFuncSpace(int o=-1) const;
virtual bool isInside(double u, double v, double w)
{
double tol = _isInsideTolerance;
if(u < (-tol) || v < (-tol) || u > ((1. + tol) - v))
return false;
return true;
}
virtual void getIntegrationPoints(int pOrder, int *npts, IntPt **pts);
virtual SPoint3 circumcenter();
private:
int edges_tri(const int edge, const int vert) const
{
static const int e[3][2] = {
{0, 1},
{1, 2},
{2, 0}
};
return e[edge][vert];
}
public:
static void registerBindings(binding *b);
};
/*
* MTriangle6
*
* 2
* |`\
* | `\
* 5 `4
* | `\
* | `\
* 0-----3----1
*
*/
class MTriangle6 : public MTriangle {
protected:
MVertex *_vs[3];
public :
MTriangle6(MVertex *v0, MVertex *v1, MVertex *v2, MVertex *v3, MVertex *v4,
MVertex *v5, int num=0, int part=0)
: MTriangle(v0, v1, v2, num, part)
{
_vs[0] = v3; _vs[1] = v4; _vs[2] = v5;
for(int i = 0; i < 3; i++) _vs[i]->setPolynomialOrder(2);
}
MTriangle6(std::vector<MVertex*> &v, int num=0, int part=0)
: MTriangle(v, num, part)
{
for(int i = 0; i < 3; i++) _vs[i] = v[3 + i];
for(int i = 0; i < 3; i++) _vs[i]->setPolynomialOrder(2);
}
~MTriangle6(){}
virtual int getPolynomialOrder() const { return 2; }
virtual int getNumVertices() const { return 6; }
virtual MVertex *getVertex(int num){ return num < 3 ? _v[num] : _vs[num - 3]; }
virtual MVertex *getVertexUNV(int num)
{
static const int map[6] = {0, 3, 1, 4, 2, 5};
return getVertex(map[num]);
}
virtual MVertex *getVertexMED(int num)
{
static const int map[6] = {0, 2, 1, 5, 4, 3};
return getVertex(map[num]);
}
virtual int getNumEdgeVertices() const { return 3; }
virtual int getNumEdgesRep();
virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(3);
MTriangle::_getEdgeVertices(num, v);
v[2] = _vs[num];
}
virtual int getNumFacesRep();
virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(6);
MTriangle::_getFaceVertices(v);
v[3] = _vs[0];
v[4] = _vs[1];
v[5] = _vs[2];
}
virtual int getTypeForMSH() const { return MSH_TRI_6; }
virtual int getTypeForUNV() const { return 92; } // thin shell parabolic triangle
//virtual int getTypeForVTK() const { return 22; }
virtual const char *getStringForPOS() const { return "ST2"; }
virtual const char *getStringForBDF() const { return "CTRIA6"; }
virtual const char *getStringForDIFF() const { return "ElmT6n2D"; }
virtual const char *getStringForINP() const { return "C2D6"; }
virtual void revert()
{
MVertex *tmp;
tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
tmp = _vs[0]; _vs[0] = _vs[2]; _vs[2] = tmp;
}
};
/*
* MTriangleN FIXME: check the plot
*
* 2
* |`\ E = order - 1;
* | `\ N = total number of vertices
* 3+2E 2+2E
* | `\ Interior vertex numbers
* ... ... for edge 0 <= i <= 2: 3+i*E to 2+(i+1)*E
* | `\ in volume : 3+3*E to N-1
* 2+3E 3+E
* | 3+3E to N-1 `\
* | `\
* 0---3--...---2+E---1
*
*/
class MTriangleN : public MTriangle {
protected:
std::vector<MVertex *> _vs;
const char _order;
public:
MTriangleN(MVertex *v0, MVertex *v1, MVertex *v2,
std::vector<MVertex*> &v, char order, int num=0, int part=0)
: MTriangle(v0, v1, v2, num, part), _vs(v), _order(order)
{
for(unsigned int i = 0; i < _vs.size(); i++) _vs[i]->setPolynomialOrder(_order);
}
MTriangleN(std::vector<MVertex*> &v, char order, int num=0, int part=0)
: MTriangle(v[0], v[1], v[2], num, part), _order(order)
{
for(unsigned int i = 3; i < v.size(); i++) _vs.push_back(v[i]);
for(unsigned int i = 0; i < _vs.size(); i++) _vs[i]->setPolynomialOrder(_order);
}
~MTriangleN(){}
virtual int getPolynomialOrder() const { return _order; }
virtual int getNumVertices() const { return 3 + _vs.size(); }
virtual MVertex *getVertex(int num){ return num < 3 ? _v[num] : _vs[num - 3]; }
virtual int getNumFaceVertices() const
{
if(_order == 3 && _vs.size() == 6) return 0;
if(_order == 3 && _vs.size() == 7) return 1;
if(_order == 4 && _vs.size() == 9) return 0;
if(_order == 4 && _vs.size() == 12) return 3;
if(_order == 5 && _vs.size() == 12) return 0;
if(_order == 5 && _vs.size() == 18) return 6;
if(_order == 6 && _vs.size() == 25) return 10;
if(_order == 7 && _vs.size() == 33) return 12;
if(_order == 8 && _vs.size() == 42) return 15;
if(_order == 9 && _vs.size() == 52) return 21;
if(_order == 10 && _vs.size() == 63) return 28;
if(_order == 6 && _vs.size() == 15) return 0;
if(_order == 7 && _vs.size() == 18) return 0;
if(_order == 8 && _vs.size() == 21) return 0;
if(_order == 9 && _vs.size() == 24) return 0;
if(_order == 10 && _vs.size() == 27) return 0;
return 0;
}
virtual int getNumEdgeVertices() const { return 3 * (_order - 1); }
virtual int getNumEdgesRep();
virtual int getNumFacesRep();
virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual void getEdgeVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(_order + 1);
MTriangle::_getEdgeVertices(num, v);
int j = 2;
const int ie = (num + 1) * (_order - 1);
for(int i = num * (_order-1); i != ie; ++i) v[j++] = _vs[i];
}
virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n);
virtual void getFaceVertices(const int num, std::vector<MVertex*> &v) const
{
v.resize(3 + _vs.size());
MTriangle::_getFaceVertices(v);
for(unsigned int i = 0; i != _vs.size(); ++i) v[i + 3] = _vs[i];
}
virtual int getTypeForMSH() const
{
if(_order == 2 && _vs.size() == 3) return MSH_TRI_6;
if(_order == 3 && _vs.size() == 6) return MSH_TRI_9;
if(_order == 3 && _vs.size() == 7) return MSH_TRI_10;
if(_order == 4 && _vs.size() == 9) return MSH_TRI_12;
if(_order == 4 && _vs.size() == 12) return MSH_TRI_15;
if(_order == 5 && _vs.size() == 12) return MSH_TRI_15I;
if(_order == 5 && _vs.size() == 18) return MSH_TRI_21;
if(_order == 6 && _vs.size() == 25) return MSH_TRI_28;
if(_order == 7 && _vs.size() == 33) return MSH_TRI_36;
if(_order == 8 && _vs.size() == 42) return MSH_TRI_45;
if(_order == 9 && _vs.size() == 52) return MSH_TRI_55;
if(_order ==10 && _vs.size() == 63) return MSH_TRI_66;
return 0;
}
virtual void revert()
{
MVertex *tmp;
tmp = _v[1]; _v[1] = _v[2]; _v[2] = tmp;
std::vector<MVertex*> inv;
inv.insert(inv.begin(), _vs.rbegin(), _vs.rend());
_vs = inv;
}
};
template <class T>
void sort3(T *t[3])
{
T *temp;
if(t[0] > t[1]){
temp = t[1];
t[1] = t[0];
t[0] = temp;
}
if(t[1] > t[2]){
temp = t[2];
t[2] = t[1];
t[1] = temp;
}
if(t[0] > t[1]){
temp = t[1];
t[1] = t[0];
t[0] = temp;
}
}
struct compareMTriangleLexicographic
{
bool operator () (MTriangle *t1, MTriangle *t2) const
{
MVertex *_v1[3] = {t1->getVertex(0), t1->getVertex(1), t1->getVertex(2)};
MVertex *_v2[3] = {t2->getVertex(0), t2->getVertex(1), t2->getVertex(2)};
sort3(_v1);
sort3(_v2);
if(_v1[0] < _v2[0]) return true;
if(_v1[0] > _v2[0]) return false;
if(_v1[1] < _v2[1]) return true;
if(_v1[1] > _v2[1]) return false;
if(_v1[2] < _v2[2]) return true;
return false;
}
};
#endif