Skip to content
Snippets Groups Projects
Commit 91e6084d authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

cleanup

parent 5758e7ea
No related branches found
No related tags found
No related merge requests found
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#include "closestPoint.h"
#include "GEntity.h"
#include "GEdge.h"
......
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _CLOSEST_POINT_H_
#define _CLOSEST_POINT_H_
#include "GmshConfig.h"
#if defined(HAVE_ANN)
#include "ANN/ANN.h"
#endif
#include "SPoint3.h"
class GEntity;
class closestPointFinder
{
......
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#include <cstdio>
#include <cmath>
#include <vector>
......@@ -7,18 +12,19 @@
#include "gmshEdge.h"
#include "Geo.h"
class discreteList {
std::vector<std::pair<SPoint3, double> > _pts;
std::vector<int> _next;
public:
int insertPoint(int pos, const SPoint3 &pt, double t) {
int insertPoint(int pos, const SPoint3 &pt, double t)
{
_pts.push_back(std::make_pair(pt, t));
_next.push_back(_next[pos + 1]);
_next[pos + 1] = _pts.size() - 1;
return _pts.size() - 1;
}
void sort(std::vector<SPoint3> &spts, std::vector<double> &ts) {
void sort(std::vector<SPoint3> &spts, std::vector<double> &ts)
{
spts.clear();
spts.reserve(_pts.size());
ts.clear();
......@@ -28,12 +34,15 @@ class discreteList {
ts.push_back(_pts[p].second);
}
}
discreteList() {
discreteList()
{
_next.push_back(-1);
}
};
static void decasteljau(double tol, discreteList &discrete, int pos, const SPoint3 &p0, const SPoint3 &p1, const SPoint3 &p2, const SPoint3 &p3, double t0, double t3)
static void decasteljau(double tol, discreteList &discrete, int pos,
const SPoint3 &p0, const SPoint3 &p1, const SPoint3 &p2,
const SPoint3 &p3, double t0, double t3)
{
SVector3 d30 = p3 - p0;
SVector3 d13 = p1 - p3;
......@@ -60,7 +69,9 @@ static void decasteljau(double tol, discreteList &discrete, int pos, const SPoin
decasteljau(tol, discrete, newpos, p0123, p123, p23, p3, t0123, t3);
}
static int discretizeBezier(double tol, discreteList &discrete, int pos, const SPoint3 pt[4], double t0, double t3, bool insertFirstPoint)
static int discretizeBezier(double tol, discreteList &discrete, int pos,
const SPoint3 pt[4], double t0, double t3,
bool insertFirstPoint)
{
if (insertFirstPoint)
pos = discrete.insertPoint(pos, pt[0], t0);
......@@ -69,7 +80,9 @@ static int discretizeBezier(double tol, discreteList &discrete, int pos, const S
return newp;
}
static int discretizeBSpline(double tol, discreteList &discrete, int pos, const SPoint3 pt[4], double t0, double t3, bool insertFirstPoint)
static int discretizeBSpline(double tol, discreteList &discrete, int pos,
const SPoint3 pt[4], double t0, double t3,
bool insertFirstPoint)
{
SPoint3 bpt[4] = {
SPoint3((pt[0] + 4 * pt[1] + pt[2]) * (1./6.)),
......@@ -80,7 +93,9 @@ static int discretizeBSpline(double tol, discreteList &discrete, int pos, const
return discretizeBezier(tol, discrete, pos, bpt, t0, t3, insertFirstPoint);
}
static int discretizeCatmullRom(double tol, discreteList &discrete, int pos, const SPoint3 pt[4], double t0, double t3, bool insertFirstPoint)
static int discretizeCatmullRom(double tol, discreteList &discrete, int pos,
const SPoint3 pt[4], double t0, double t3,
bool insertFirstPoint)
{
SPoint3 bpt[4] = {
pt[1],
......@@ -98,7 +113,8 @@ static inline SPoint3 curveGetPoint(Curve *c, int i)
return SPoint3(v->Pos.X, v->Pos.Y, v->Pos.Z);
}
static void discretizeCurve(Curve *c, double tol, std::vector<SPoint3> &pts, std::vector<double> &ts)
static void discretizeCurve(Curve *c, double tol, std::vector<SPoint3> &pts,
std::vector<double> &ts)
{
discreteList discrete;
switch(c->Typ) {
......@@ -190,7 +206,8 @@ static void discretizeCurve(Curve *c, double tol, std::vector<SPoint3> &pts, std
discrete.sort(pts, ts);
}
void gmshEdge::discretize(double tol, std::vector<SPoint3> &dpts, std::vector<double> &ts)
void gmshEdge::discretize(double tol, std::vector<SPoint3>
&dpts, std::vector<double> &ts)
{
discretizeCurve(c, tol, dpts, ts);
}
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#include "BGMBase.h"
#include <iostream>
#include "OS.h"
#include "GPoint.h"
#include "GFace.h"
#include "GmshDefines.h"
#include "MElementOctree.h"
void BGMBase::export_scalar(const std::string &filename, const DoubleStorageType &_whatToPrint) const
void BGMBase::export_scalar(const std::string &filename,
const DoubleStorageType &_whatToPrint) const
{
FILE *f = Fopen (filename.c_str(),"w");
fprintf(f,"View \"Background Mesh\"{\n");
......@@ -56,7 +61,8 @@ void BGMBase::export_scalar(const std::string &filename, const DoubleStorageType
fclose(f);
}
void BGMBase::export_vector(const std::string &filename, const VectorStorageType &_whatToPrint) const
void BGMBase::export_vector(const std::string &filename,
const VectorStorageType &_whatToPrint) const
{
FILE *f = Fopen (filename.c_str(),"w");
fprintf(f,"View \"Background Mesh\"{\n");
......@@ -108,7 +114,8 @@ void BGMBase::export_vector(const std::string &filename, const VectorStorageType
}
void BGMBase::export_tensor_as_vectors(const std::string &filename, const TensorStorageType &_whatToPrint)const
void BGMBase::export_tensor_as_vectors(const std::string &filename,
const TensorStorageType &_whatToPrint) const
{
FILE *f = Fopen (filename.c_str(),"w");
fprintf(f,"View \"Background Mesh\"{\n");
......@@ -120,20 +127,22 @@ void BGMBase::export_tensor_as_vectors(const std::string &filename, const Tensor
v = it->first;
GPoint p = get_GPoint_from_MVertex(v);
for (int i=0;i<3;i++){
fprintf(f,"%s(%g,%g,%g){%g,%g,%g};\n",s,p.x(),p.y(),p.z(),(it->second)(0,i),(it->second)(1,i),(it->second)(2,i));
fprintf(f,"%s(%g,%g,%g){%g,%g,%g};\n",s,p.x(),p.y(),p.z(),-(it->second)(0,i),-(it->second)(1,i),-(it->second)(2,i));
fprintf(f,"%s(%g,%g,%g){%g,%g,%g};\n",s,p.x(),p.y(),p.z(),
(it->second)(0,i),(it->second)(1,i),(it->second)(2,i));
fprintf(f,"%s(%g,%g,%g){%g,%g,%g};\n",s,p.x(),p.y(),p.z(),
-(it->second)(0,i),-(it->second)(1,i),-(it->second)(2,i));
}
}
fprintf(f,"};\n");
fclose(f);
}
BGMBase::BGMBase(int dim,GEntity *_gf):octree(NULL),gf(_gf), DIM(dim), order(1)
{
}
BGMBase::~BGMBase(){
BGMBase::~BGMBase()
{
}
bool BGMBase::inDomain (double u, double v, double w)
......@@ -146,7 +155,8 @@ const MElement* BGMBase::findElement(double u, double v, double w, bool strict)
return (getOctree()->find(u, v, w, DIM, strict));
}
std::vector<double> BGMBase::get_field_value(double u, double v, double w, const VectorStorageType &data)
std::vector<double> BGMBase::get_field_value(double u, double v, double w,
const VectorStorageType &data)
{
MElement *e = const_cast<MElement*>(findElement(u, v, w ));
if (!e) return std::vector<double>(3,-1000.);
......@@ -157,12 +167,14 @@ std::vector<double> BGMBase::get_field_value(double u, double v, double w, const
for (int j=0;j<3;j++){
std::vector<double> values(e->getNumVertices());
for (int i=0;i<e->getNumVertices();i++) values[i]=val[i][j];
res[j] = e->interpolate(&values[0], element_uvw[0], element_uvw[1], element_uvw[2], 1, order);
res[j] = e->interpolate(&values[0], element_uvw[0], element_uvw[1],
element_uvw[2], 1, order);
}
return res;
}
double BGMBase::get_field_value(double u, double v, double w, const DoubleStorageType &data)
double BGMBase::get_field_value(double u, double v, double w,
const DoubleStorageType &data)
{
MElement *e = const_cast<MElement*>(findElement(u, v, w));
if (!e) return -1000.;
......@@ -172,7 +184,8 @@ double BGMBase::get_field_value(double u, double v, double w, const DoubleStorag
for (int i=0;i<e->getNumVertices();i++)
values[i]=val[i];
return e->interpolate(&values[0], element_uvw[0], element_uvw[1], element_uvw[2], 1, order);
return e->interpolate(&values[0], element_uvw[0], element_uvw[1],
element_uvw[2], 1, order);
}
double BGMBase::size(double u, double v, double w)
......@@ -185,7 +198,8 @@ double BGMBase::size(const MVertex *v)
return get_nodal_value(v,sizeField);
}
std::vector<double> BGMBase::get_nodal_value(const MVertex *v,const VectorStorageType &data)const
std::vector<double> BGMBase::get_nodal_value(const MVertex *v,
const VectorStorageType &data) const
{
VectorStorageType::const_iterator itfind = data.find(const_cast<MVertex*>(v));
if (itfind==data.end()){
......@@ -205,7 +219,8 @@ double BGMBase::get_nodal_value(const MVertex *v,const DoubleStorageType &data)c
return itfind->second;
}
std::vector<std::vector<double> > BGMBase::get_nodal_values(const MElement *e,const VectorStorageType &data)const
std::vector<std::vector<double> >
BGMBase::get_nodal_values(const MElement *e,const VectorStorageType &data) const
{
std::vector<std::vector<double> > res(e->getNumVertices());
......@@ -217,7 +232,8 @@ std::vector<std::vector<double> > BGMBase::get_nodal_values(const MElement *e,co
return res;
}
std::vector<double> BGMBase::get_nodal_values(const MElement *e,const DoubleStorageType &data)const
std::vector<double> BGMBase::get_nodal_values(const MElement *e,
const DoubleStorageType &data) const
{
std::vector<double> res(e->getNumVertices(),0.);
......@@ -226,7 +242,8 @@ std::vector<double> BGMBase::get_nodal_values(const MElement *e,const DoubleStor
return res;
}
std::vector<double> BGMBase::get_element_uvw_from_xyz (const MElement *e, double x, double y,double z) const
std::vector<double> BGMBase::get_element_uvw_from_xyz (const MElement *e, double x,
double y, double z) const
{
double element_uvw[3];
double xyz[3] = {x, y, z};
......@@ -255,7 +272,3 @@ GEntity* BGMBase::getBackgroundGEntity()
{
return gf;
}
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _BACKGROUND_MESH_TOOLS_H_
#define _BACKGROUND_MESH_TOOLS_H_
......@@ -10,24 +13,18 @@ class GVertex;
class GEdge;
class GEntity;
SMetric3 buildMetricTangentToCurve(SVector3 &t, double l_t, double l_n);
SMetric3 buildMetricTangentToSurface (SVector3 &t1, SVector3 &t2, double l_t1, double l_t2, double l_n);
SMetric3 buildMetricTangentToSurface(SVector3 &t1, SVector3 &t2, double l_t1,
double l_t2, double l_n);
double BGM_MeshSize(GEntity *ge, double U, double V, double X, double Y, double Z);
SMetric3 BGM_MeshMetric(GEntity *ge, double U, double V, double X, double Y, double Z);
bool Extend1dMeshIn2dSurfaces();
bool Extend2dMeshIn3dVolumes();
SMetric3 max_edge_curvature_metric(const GVertex *gv);
SMetric3 max_edge_curvature_metric(const GEdge *ge, double u, double &l);
SMetric3 metric_based_on_surface_curvature(const GFace *gf, double u, double v, bool surface_isotropic = false,double d_normal = 1.e12,double d_tangent_max = 1.e12);
//
//static double LC_MVertex_CURV(GEntity *ge, double U, double V);
//SMetric3 LC_MVertex_CURV_ANISO(GEntity *ge, double U, double V);
//static double LC_MVertex_PNTS(GEntity *ge, double U, double V);
//static double max_edge_curvature(const GVertex *gv);
//static double max_surf_curvature(const GEdge *ge, double u);
//static SMetric3 metric_based_on_surface_curvature(const GEdge *ge, double u, bool iso_surf);
//static SMetric3 metric_based_on_surface_curvature(const GVertex *gv, bool iso_surf);
SMetric3 metric_based_on_surface_curvature(const GFace *gf, double u, double v,
bool surface_isotropic = false,
double d_normal = 1.e12,
double d_tangent_max = 1.e12);
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _FIELD_PYTHON_H_
#define _FIELD_PYTHON_H_
#include "Field.h"
#include "Python.h"
class FieldPython : public Field
{
PyObject *_callback;
......@@ -56,4 +62,5 @@ class FieldPython : public Field
}
}
};
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _FILTER_OVERLAPPING_ELEMENTS_
#define _FILTER_OVERLAPPING_ELEMENTS_
#include <map>
#include <vector>
class MElement;
class MPrism;
class MHexahedron;
class MTriangle;
class MQuadrangle;
void filterOverlappingElements(std::vector<MElement*> &,
std::map<MElement*,std::vector<MElement*> > &_elemColumns,
std::map<MElement*,MElement*> &_toFirst);
......@@ -20,4 +28,5 @@ void filterOverlappingElements (std::vector<MTriangle*> &blTris,
std::map<MElement*,MElement*> &_toFirst);
void filterColumns(std::vector<MElement*> &elem,
std::map<MElement*,std::vector<MElement*> > &_elemColumns);
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#include "approximationError.h"
#include "MElement.h"
double approximationError (simpleFunction<double> &f, MElement *e)
{
std::vector<double> VALS(e->getNumVertices());
......
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _APPROXIMATION_ERROR_
#define _APPROXIMATION_ERROR_
#include <map>
#include "simpleFunction.h"
class MElement;
// computes E such as
// E^2 = \int_e [C_e(f) - f]^2 de
// where C_e(f) is clement's interpolation operator of f on e
double approximationError(simpleFunction<double> &f, MElement *e);
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _DISCRETE_FRECHET_DISTANCE_
#define _DISCRETE_FRECHET_DISTANCE_
#include <vector>
#include "SPoint3.h"
double discreteFrechetDistance (const std::vector<SPoint3> &P,
const std::vector<SPoint3> &Q);
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
/*
compute the hausdorff distance between two polygonal curves
in n*m time where n and m are the nb of points of the
polygonal curves
*/
#include "SVector3.h"
#include "hausdorffDistance.h"
......@@ -20,7 +26,8 @@ static double intersect (SPoint3 &q, SVector3 &n, // a plane
return t;
}
static double projection (SPoint3 &p1, SPoint3 &p2, SPoint3 &q, SPoint3 &result){
static double projection (SPoint3 &p1, SPoint3 &p2, SPoint3 &q, SPoint3 &result)
{
// x = p1 + t (p2 - p1)
// (x - q) * (p2 - p1) = 0
// (p1 + t (p2 - p1) - q) (p2 - p1) = 0
......@@ -31,7 +38,8 @@ static double projection (SPoint3 &p1, SPoint3 &p2, SPoint3 &q, SPoint3 &result)
return t;
}
static SPoint3 closestPoint (SPoint3 &p1, SPoint3 &p2, SPoint3 &q){
static SPoint3 closestPoint (SPoint3 &p1, SPoint3 &p2, SPoint3 &q)
{
double result;
const double t = projection (p1,p2,q,result);
if (t >= 0.0 && t <= 1.0)return result;
......@@ -39,8 +47,8 @@ static SPoint3 closestPoint (SPoint3 &p1, SPoint3 &p2, SPoint3 &q){
return p2;
}
double closestPoint (const std::vector<SPoint3> &P, const SPoint3 &p, SPoint3 & result){
double closestPoint (const std::vector<SPoint3> &P, const SPoint3 &p, SPoint3 & result)
{
double closestDistance = 1.e22;
for (unsigned int i=1;i<P.size();i++){
SPoint3 q = closestPoint (P[i-1],P[i],p);
......@@ -57,7 +65,8 @@ double closestPoint (const std::vector<SPoint3> &P, const SPoint3 &p, SPoint3 &
// of angle bissectors of Q with P
double oneSidedHausdorffDistance (const std::vector<SPoint3> &P,
const std::vector<SPoint3> &Q,
SPoint3 &p1, SPoint3 &p2){
SPoint3 &p1, SPoint3 &p2)
{
const double hausdorffDistance = 0.0;
// first test the points
......@@ -109,7 +118,8 @@ double oneSidedHausdorffDistance (const std::vector<SPoint3> &P,
}
double hausdorffDistance (const std::vector<SPoint3> &P,
const std::vector<SPoint3> &Q){
const std::vector<SPoint3> &Q)
{
return std::max(oneSidedHausdorffDistance (P,Q),
oneSidedHausdorffDistance (Q,P));
}
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#include "miniBasis.h"
#include "BasisFactory.h"
miniBasisTri::miniBasisTri()
{
type = MSH_TRI_MINI;
......
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _MINI_BASIS_H_
#define _MINI_BASIS_H_
#include "polynomialBasis.h"
// mini is NOT a real nodal basis but in GMSH, only the nodal basis have closure and mini have closure so...
// mini is NOT a real nodal basis but in GMSH, only the nodal basis have closure
// and mini have closure so...
class miniBasisTri : public polynomialBasis {
public:
miniBasisTri();
......@@ -10,4 +18,5 @@ class miniBasisTet : public polynomialBasis {
public:
miniBasisTet();
};
#endif
// Gmsh - Copyright (C) 1997-2015 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
#ifndef _SIMPLE_FUNCTION_PYTHON_H_
#define _SIMPLE_FUNCTION_PYTHON_H_
#include "Python.h"
#include "simpleFunction.h"
......@@ -34,4 +40,5 @@ class simpleFunctionPython : public simpleFunction<double> {
return r;
}
};
#endif
......@@ -186,7 +186,8 @@ information.
<h2><a name="Screenshots"></a>Screenshots</h2>
<a href="gallery/screenshot.png"><img src="gallery/thumbnail.png"
alt="Screenshot thumbnail"></a>
alt="Screenshot thumbnail"></a>&nbsp;<a href="gallery/screenshot2.png"><img
src="gallery/thumbnail2.png" alt="Screenshot thumbnail"></a>
<ul>
<li>Sample STEP/BREP models:
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment