Skip to content
Snippets Groups Projects
Commit 0f09ea5b authored by Paul-Emile Bernard's avatar Paul-Emile Bernard
Browse files

Adding technique 5 to meshMetric: linear interpolation of h in band around...

Adding technique 5 to meshMetric: linear interpolation of h in band around iso, instead of linear interpolation of the metric (1/h^2)
parent 9cada17f
No related branches found
No related tags found
No related merge requests found
......@@ -548,7 +548,7 @@ int GModel::adaptMesh(std::vector<int> technique, std::vector<simpleFunction<dou
fields->reset();
meshMetric *metric = new meshMetric(this);
for (int imetric=0;imetric<technique.size();imetric++){;
for (int imetric=0;imetric<technique.size();imetric++){
metric->addMetric(technique[imetric], f[imetric], parameters[imetric]);
}
fields->setBackgroundField(metric);
......@@ -561,9 +561,9 @@ int GModel::adaptMesh(std::vector<int> technique, std::vector<simpleFunction<dou
opt_mesh_algo3d(0, GMSH_SET, 7.0); //mmg3d
opt_mesh_lc_from_points(0, GMSH_SET, 0.0); //do not mesh lines with lc
std::for_each(firstEdge(), lastEdge(), deMeshGEdge());
std::for_each(firstFace(), lastFace(), deMeshGFace());
std::for_each(firstRegion(), lastRegion(), deMeshGRegion());
std::for_each(firstFace(), lastFace(), deMeshGFace());
std::for_each(firstEdge(), lastEdge(), deMeshGEdge());
GenerateMesh(this, getDim());
nbElems = getNumMeshElements();
......@@ -606,7 +606,7 @@ int GModel::adaptMesh(std::vector<int> technique, std::vector<simpleFunction<dou
fields->reset();
meshMetric *metric = new meshMetric(this);
for (int imetric=0;imetric<technique.size();imetric++){;
for (int imetric=0;imetric<technique.size();imetric++){
metric->addMetric(technique[imetric], f[imetric], parameters[imetric]);
}
fields->setBackgroundField(metric);
......
......@@ -396,6 +396,7 @@ class GModel
// parameters[0] = thickness of the interface in the positive ls direction (mandatory)
// parameters[4] = thickness of the interface in the negative ls direction (=parameters[0] if not specified)
// parameters[3] = the required minimum number of elements to represent a circle - used for curvature-based metric (default: = 15)
// 5) Same as 4, except that the transition in band E uses linear interpolation of h, instead of linear interpolation of metric
// The algorithm first generate a mesh if no one is available
// In this first attempt, only the highest dimensional mesh is adapted, which is ok if
......
......@@ -73,6 +73,7 @@ void meshMetric::addMetric(int technique, simpleFunction<double> *fct, std::vect
_E = parameters[0];
_E_moins = (parameters.size() >= 5) ? parameters[4] : -parameters[0];
if (_E_moins>0.) _E_moins *= -1.;
_epsilon = parameters[0];
_technique = (MetricComputationTechnique)technique;
_Np = (parameters.size() >= 4) ? parameters[3] : 15.;
......@@ -342,7 +343,7 @@ void meshMetric::computeMetric(){
// else H = hfrey;
H = hfrey;
}
else if (_technique == meshMetric::EIGENDIRECTIONS ){
else if ((_technique == meshMetric::EIGENDIRECTIONS )||(_technique == meshMetric::EIGENDIRECTIONS_LINEARINTERP_H )){
double metric_value_hmax = 1./hmax/hmax;
SVector3 gr = grads[ver];
double norm = gr.normalize();
......@@ -369,12 +370,16 @@ void meshMetric::computeMetric(){
eigenvals_curvature.push_back(fabs(lambda_hessian(2))*un_sur_epsilon);
// metric distance-based eigenvalue, in gradient direction
double metric_value_hmin = 1./hmin/hmin;
// 1/ linear interpolation between h_min and h_max...
//double h_dist = hmin + ((hmax-hmin)/_E)*dist;// the charcteristic element size in the normal direction - linear interp between hmin and hmax
//eigenval_direction = 1./h_dist/h_dist;
// 2/ ... or linear interpolation between 1/h_min^2 and 1/h_max^2
double eigenval_direction;
if (_technique==meshMetric::EIGENDIRECTIONS_LINEARINTERP_H){
double h_dist = hmin + ((hmax-hmin)/_E)*dist;// the charcteristic element size in the normal direction - linear interp between hmin and hmax
eigenval_direction = 1./h_dist/h_dist;
}
else if(_technique==meshMetric::EIGENDIRECTIONS){
// ... or linear interpolation between 1/h_min^2 and 1/h_max^2
double maximum_distance = (signed_dist>0.) ? _E : fabs(_E_moins);
double eigenval_direction = metric_value_hmin + ((metric_value_hmax-metric_value_hmin)/maximum_distance)*dist;
eigenval_direction = metric_value_hmin + ((metric_value_hmax-metric_value_hmin)/maximum_distance)*dist;
}
// now, consider only three eigenvectors (i.e. (grad(LS), u1, u2) with u1 and u2 orthogonal vectors in the tangent plane to the LS)
// and imposing directly eigenvalues and directions, instead of intersecting the two metrics based on direction and curvature
......@@ -413,7 +418,7 @@ void meshMetric::computeMetric(){
it++;
}
if (_technique != meshMetric::EIGENDIRECTIONS ){
if (_technique != meshMetric::EIGENDIRECTIONS && _technique!=meshMetric::EIGENDIRECTIONS_LINEARINTERP_H){
//See paper Hachem and Coupez:
//Finite element solution to handle complex heat and fluid flows in
......
......@@ -14,7 +14,7 @@ class STensor3;
/**Anisotropic mesh size field based on a metric */
class meshMetric: public Field {
public:
typedef enum {LEVELSET=1,HESSIAN=2, FREY=3, EIGENDIRECTIONS=4} MetricComputationTechnique;
typedef enum {LEVELSET=1,HESSIAN=2, FREY=3, EIGENDIRECTIONS=4, EIGENDIRECTIONS_LINEARINTERP_H=5} MetricComputationTechnique;
private:
// intersect all metrics added in "setOfMetrics", preserve eigendirections of the "most anisotropic" metric
void intersectMetrics();
......@@ -63,6 +63,7 @@ class meshMetric: public Field {
// parameters[0] = thickness of the interface in the positive ls direction (mandatory)
// parameters[4] = thickness of the interface in the negative ls direction (default: =parameters[0] if not specified)
// parameters[3] = the required minimum number of elements to represent a circle - used for curvature-based metric (default: = 15)
// 5: same as 4, except that the transition in band E uses linear interpolation of h, instead of linear interpolation of metric
void addMetric(int technique, simpleFunction<double> *fct, std::vector<double> parameters);
inline SMetric3 metricAtVertex (MVertex* v) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment