diff --git a/Graphics/Geom.cpp b/Graphics/Geom.cpp index 719e5990bc6c685076711f72ba8d7cf3c5a3c531..30c5cffdc8010a7caf9097a2232e54b15116ac52 100644 --- a/Graphics/Geom.cpp +++ b/Graphics/Geom.cpp @@ -1,4 +1,4 @@ -// $Id: Geom.cpp,v 1.66 2004-06-22 17:34:10 geuzaine Exp $ +// $Id: Geom.cpp,v 1.67 2004-06-23 03:57:43 geuzaine Exp $ // // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle // @@ -285,8 +285,6 @@ void Draw_Triangulated_Surface(Surface * s) } } -void Get_SurfaceNormal(Surface *s, double n[3]); - void Draw_Plane_Surface(Surface * s) { int i, j, k; diff --git a/Mesh/2D_Mesh.cpp b/Mesh/2D_Mesh.cpp index 0ddb2885287734e6c8fcb0e379a1073e75b883e3..e42f166d3ba4e87e497ee1901c93c803a130d0a1 100644 --- a/Mesh/2D_Mesh.cpp +++ b/Mesh/2D_Mesh.cpp @@ -1,4 +1,4 @@ -// $Id: 2D_Mesh.cpp,v 1.63 2004-06-22 22:10:11 geuzaine Exp $ +// $Id: 2D_Mesh.cpp,v 1.64 2004-06-23 03:57:43 geuzaine Exp $ // // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle // @@ -773,23 +773,21 @@ void ActionEndTheCurve(void *a, void *b) End_Curve(c); } -void ActionInvertTriQua(void *a, void *b) +void ActionInvertTri(void *a, void *b) { Simplex *s = *(Simplex **) a; - Vertex *tmp; - if(s->V[3]){ - tmp = s->V[1]; - s->V[1] = s->V[3]; - s->V[3] = tmp; - } - else if(s->V[2]){ - tmp = s->V[1]; - s->V[1] = s->V[2]; - s->V[2] = tmp; - } + Vertex *tmp = s->V[1]; + s->V[1] = s->V[2]; + s->V[2] = tmp; } -int isPointOnPlanarSurface(Surface * S, double X, double Y, double Z, double n[3]); +void ActionInvertQua(void *a, void *b) +{ + Quadrangle *q = *(Quadrangle **) a; + Vertex *tmp = q->V[1]; + q->V[1] = q->V[3]; + q->V[3] = tmp; +} int isMiddlePointOnPlanarSurface(Surface *s, Vertex *v1, Vertex *v2) { @@ -858,15 +856,19 @@ void Get_SurfaceNormal(Surface *s, double n[3]) } } +// Horrible (tm) hack to orient the elements correctly. This is +// *definitely* not the best way to do it, but I don't have time to +// look into the issue right now. + void ReOrientSurfaceMesh(Surface *s) { - // Horrible (tm) hack to orient the elements correctly. This is - // *definitely* not the best way to do it, but I don't have time - // to look into the issue right now. + double t1[3], t2[3], n1[3], n2[3], res; Simplex *simp; + Quadrangle *quad; + + Get_SurfaceNormal(s, n1); + if(Tree_Right(s->Simplexes, &simp)){ - double t1[3], t2[3], n1[3], n2[3], res; - Get_SurfaceNormal(s, n1); t1[0] = simp->V[1]->Pos.X - simp->V[0]->Pos.X; t1[1] = simp->V[1]->Pos.Y - simp->V[0]->Pos.Y; t1[2] = simp->V[1]->Pos.Z - simp->V[0]->Pos.Z; @@ -875,19 +877,28 @@ void ReOrientSurfaceMesh(Surface *s) t2[2] = simp->V[2]->Pos.Z - simp->V[0]->Pos.Z; prodve(t1, t2, n2); norme(n2); - /* - printf("n1=%g %g %g\n", n1[0], n1[1], n1[2]); - printf("n2=%g %g %g (elt: (%g,%g,%g) (%g,%g,%g) (%g,%g,%g))\n", - n2[0], n2[1], n2[2], - simp->V[0]->Pos.X, simp->V[0]->Pos.Y, simp->V[0]->Pos.Z, - simp->V[1]->Pos.X, simp->V[1]->Pos.Y, simp->V[1]->Pos.Z, - simp->V[2]->Pos.X, simp->V[2]->Pos.Y, simp->V[2]->Pos.Z); - */ prosca(n1, n2, &res); if(res < 0.0){ - Msg(DEBUG, "Inverting orientation of elements in %s surface %d (res = %g)", + Msg(DEBUG, "Inverting triangles in %s surface %d (res = %g)", (s->Typ == MSH_SURF_PLAN) ? "Plane" : "NonPlane", s->Num, res); - Tree_Action(s->Simplexes, ActionInvertTriQua); + Tree_Action(s->Simplexes, ActionInvertTri); + } + } + + if(Tree_Right(s->Quadrangles, &quad)){ + t1[0] = quad->V[1]->Pos.X - quad->V[0]->Pos.X; + t1[1] = quad->V[1]->Pos.Y - quad->V[0]->Pos.Y; + t1[2] = quad->V[1]->Pos.Z - quad->V[0]->Pos.Z; + t2[0] = quad->V[2]->Pos.X - quad->V[0]->Pos.X; + t2[1] = quad->V[2]->Pos.Y - quad->V[0]->Pos.Y; + t2[2] = quad->V[2]->Pos.Z - quad->V[0]->Pos.Z; + prodve(t1, t2, n2); + norme(n2); + prosca(n1, n2, &res); + if(res < 0.0){ + Msg(DEBUG, "Inverting quads in %s surface %d (res = %g)", + (s->Typ == MSH_SURF_PLAN) ? "Plane" : "NonPlane", s->Num, res); + Tree_Action(s->Quadrangles, ActionInvertQua); } } } @@ -929,18 +940,12 @@ void Maillage_Surface(void *data, void *dum) if(MeshTransfiniteSurface(s) || MeshEllipticSurface(s) || MeshCylindricalSurface(s) || - MeshParametricSurface(s)) { + MeshParametricSurface(s) || + Extrude_Mesh(s)) { Tree_Action(THEM->Points, PutVertex_OnSurf); Tree_Action(s->Vertices, PutVertex_OnSurf); Tree_Action(s->Vertices, Add_In_Mesh); } - else if(Extrude_Mesh(s)) { - Tree_Action(THEM->Points, PutVertex_OnSurf); - Tree_Action(s->Vertices, PutVertex_OnSurf); - Tree_Action(s->Vertices, Add_In_Mesh); - // FIXME: big hack - //ReOrientSurfaceMesh(s); - } else{ int TypSurface = s->Typ; Plan_Moyen(pS, dum); diff --git a/Mesh/3D_Extrude.cpp b/Mesh/3D_Extrude.cpp index 6920dd68ec9311739b321d779a89f7e1494cf696..f98b655dd75d2982c3dcc246165cbf8b46f2a045 100644 --- a/Mesh/3D_Extrude.cpp +++ b/Mesh/3D_Extrude.cpp @@ -1,4 +1,4 @@ -// $Id: 3D_Extrude.cpp,v 1.80 2004-06-08 02:10:32 geuzaine Exp $ +// $Id: 3D_Extrude.cpp,v 1.81 2004-06-23 03:57:43 geuzaine Exp $ // // Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle // @@ -1046,6 +1046,9 @@ int Extrude_Mesh(Surface * s) Tree_Action(s->Simplexes, AddSimVertsInSurf); Tree_Action(s->Quadrangles, AddQuadVertsInSurf); + // FIXME: big hack + ReOrientSurfaceMesh(s); + return true; } diff --git a/Mesh/Mesh.h b/Mesh/Mesh.h index 4f3cb679be8c4458c248ec019ab8bb0a22684b96..181cc9c19eeecf2dd18c8b15e53bca20c2502ffa 100644 --- a/Mesh/Mesh.h +++ b/Mesh/Mesh.h @@ -476,6 +476,9 @@ void Projette_Plan_Moyen(void *a, void *b); void Projette_Inverse(void *a, void *b); void Freeze_Vertex(void *a, void *b); void deFreeze_Vertex(void *a, void *b); +int isPointOnPlanarSurface(Surface * S, double X, double Y, double Z, double n[3]); +void Get_SurfaceNormal(Surface *s, double n[3]); +void ReOrientSurfaceMesh(Surface *s); double Lc_XYZ(double X, double Y, double Z, Mesh * m); void ActionLiss(void *data, void *dummy); diff --git a/TODO b/TODO index 0c7cfb233c177fc4a31321dfcaae5dc0b7042aa0..a25a41f5ccd5230e6067d2b8e7d9865d75c67899 100644 --- a/TODO +++ b/TODO @@ -1,4 +1,4 @@ -$Id: TODO,v 1.49 2004-06-02 02:22:13 geuzaine Exp $ +$Id: TODO,v 1.50 2004-06-23 03:57:43 geuzaine Exp $ add ternary operator and <,>,<=,>=,== tests in MathEval @@ -34,6 +34,12 @@ interactive use?) ******************************************************************** +Find a better solution for the orientation of the surface mesh +elements. The current "a posteriori re-orientation" solution is a +real crime against good taste :-) + +******************************************************************** + Yves Krahenbuhl wrote: > Lors de la creation des elements du 2eme ordre, et selon la courbure @@ -45,7 +51,7 @@ Yves Krahenbuhl wrote: ******************************************************************** -The "Symmetry" operation should be renamed "Reflection" +The "Symmetry" operation should be renamed "Reflection" (?) ******************************************************************** @@ -55,7 +61,8 @@ Attractors in the 2D aniso algo are extremely buggy Memory leaks, memory leaks -- start with mesh_domain() and the parser +- start with mesh_domain() and the parser (update: parser should be + mostly OK now) - check all calls to Tree_Replace: we shouldn't use it with trees of pointers, since we loose the original pointer when we actually do a @@ -90,6 +97,8 @@ etc. Include Tetgen (the same way we did it with Triangle--the license is the same). +******************************************************************** + Include the 3D frontal algo from Netgen? (it's LGPL now...) ******************************************************************** @@ -162,68 +171,3 @@ Remarque finale: le calcul des valeurs et vecteurs propres est à partir des v.p.), et cela nous rappelle cette belle expression "No free meal" ! -******************************************************************** - -Stockage et acces aux normales en post-pro: generer les noeuds pour -acces plus rapide? - -Idem pour le maillage. - -******************************************************************** - -Orientation des surfaces - -mesh: il faut vraiment que je regarde la relation entre l'orientation -des surfaces geometriques et les l'orientation des elements du -maillage... - -Subject: Re: [Gmsh] surface normals: some point inwards, but most -point in the negative direction - Date: Thu, 21 Nov 2002 05:59:56 -0500 - From: "Matthias G. Imhof" <mgi@vt.edu> - Organization: Virginia Tech: Seismic Reservoir Characterization -Laboratory - To: Christophe Geuzaine <geuzaine@acm.caltech.edu> - References: 1 , 2 - -Christophe, - -Thanks for a nice tool! I finally extended my boundary-element method -to 3D and find your mesh generator perfect to make sure the boundary -points are placed uniformly on a sphere instead of clustering up at -the poles. - -My only gripe is with the visualization of surface normals in -gmsh. While Physical Surface allows change of the normal directions in -the mesh file, the normals are not changed in the visualization which -makes identification of wrong polarity surfaces difficult. - -Matthias - -Christophe Geuzaine wrote: -> Matthias G. Imhof wrote: -> ->> I defined a simple rectangular block and noticed that normals in the ->> x-direction actually change signs depending on the surface to point ->> into the block, while the other normals point in the negative y- or ->> z-directions pointing into the block for some surfaces and pointing ->> outwards for others. ->> ->> How can I change this model to have all normals point out of the block? -> -> -> The only way is to define Physical Surfaces with the appropriate sign, -> e.g. "Physical Surface(1) = {23,1,20,3,2,-16};". This will reverse the -> orientation of the surface 16 in the output mesh file. -> ->> ->> A second question for the same problem. If I generate a quantity by ->> duplicate, how can I lates access this quantity? I tried to assign it ->> to, e.g., a plane surface but that did not work. -> -> -> There is unfortunately no way to access that information at the moment -> without using the GUI... -> -> Christophe ->