Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
#include <math.h>
#include "MElement.h"
#include "GEntity.h"
int MElement::_globalNum = 0;
static double dist(MVertex *v1, MVertex *v2)
{
double dx = v1->x() - v2->x();
double dy = v1->y() - v2->y();
double dz = v1->z() - v2->z();
return sqrt(dx * dx + dy * dy + dz * dz);
}
double MElement::minEdge()
{
double m = 1.e25;
MVertex *v[2];
for(int i = 0; i < getNumEdges(); i++){
getEdge(i, v);
m = std::min(m, dist(v[0], v[1]));
}
return m;
}
double MElement::maxEdge()
{
double m = 0.;
MVertex *v[2];
for(int i = 0; i < getNumEdges(); i++){
getEdge(i, v);
m = std::max(m, dist(v[0], v[1]));
}
return m;
}
double MElement::rhoShapeMeasure()
{
double min = minEdge();
double max = maxEdge();
if(max)
return min / max;
else
return 0.;
}
void MElement::cog(double &x, double &y, double &z)
{
x = y = z = 0.;
int n = getNumVertices();
for(int i = 0; i < n; i++) {
MVertex *v = getVertex(i);
x += v->x();
y += v->y();
z += v->z();
}
x /= (double)n;
y /= (double)n;
z /= (double)n;
}
void MElement::writeMSH(FILE *fp, double version, int num, int elementary,
int physical)
{
int n = getNumVertices();
fprintf(fp, "%d %d", num ? num : _num, getTypeForMSH());
if(version < 2.0)
fprintf(fp, " %d %d %d", physical, elementary, n);
else
fprintf(fp, " 3 %d %d %d", physical, elementary, _partition);
for(int i = 0; i < n; i++){
fprintf(fp, " %d", getVertex(i)->getNum());
}
fprintf(fp, " \n");
}
void MElement::writePOS(FILE *fp)
{
int n = getNumVertices();
double gamma = gammaShapeMeasure();
double eta = etaShapeMeasure();
double rho = rhoShapeMeasure();
fprintf(fp, "%s(", getStringForPOS());
for(int i = 0; i < n; i++){
if(i) fprintf(fp, ",");
fprintf(fp, "%g,%g,%g", getVertex(i)->x(), getVertex(i)->y(), getVertex(i)->z());
}
fprintf(fp, "){");
for(int i = 0; i < n; i++)
fprintf(fp, "%d,", getVertex(i)->onWhat()->tag());
for(int i = 0; i < n; i++)
fprintf(fp, "%d,", getNum());
for(int i = 0; i < n; i++)
fprintf(fp, "%g,", gamma);
for(int i = 0; i < n; i++)
fprintf(fp, "%g,", eta);
for(int i = 0; i < n; i++){
if(i == n - 1)
fprintf(fp, "%g", rho);
else
fprintf(fp, "%g,", rho);
}
fprintf(fp, "};\n");
}