-// Copyright (C) 2007-2008 CEA/DEN, EDF R&D
+// Copyright (C) 2007-2011 CEA/DEN, EDF R&D
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
//
// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
+
// ---
// File : BLSURFPlugin_BLSURF.cxx
// Authors : Francis KLOSS (OCC) & Patrick LAUG (INRIA) & Lioka RAZAFINDRAZAKA (CEA)
#include "BLSURFPlugin_Hypothesis.hxx"
extern "C"{
-#include "distene/blsurf.h"
#include <distene/api.h>
+#include <distene/blsurf.h>
}
#include <structmember.h>
+#include <Basics_Utils.hxx>
+
#include <SMESH_Gen.hxx>
#include <SMESH_Mesh.hxx>
#include <SMESH_ControlsDef.hxx>
+#include <SMDSAbs_ElementType.hxx>
+#include "SMESHDS_Group.hxx"
+#include "SMESH_Group.hxx"
#include <utilities.h>
#include <set>
#include <cstdlib>
+// OPENCASCADE includes
#include <BRep_Tool.hxx>
#include <TopExp.hxx>
#include <TopExp_Explorer.hxx>
+#include <TopoDS.hxx>
#include <NCollection_Map.hxx>
-#include <Standard_ErrorHandler.hxx>
-
-extern "C"{
-#include "distene/blsurf.h"
-#include <distene/api.h>
-}
#include <Geom_Surface.hxx>
#include <Handle_Geom_Surface.hxx>
#include <Handle_Geom2d_Curve.hxx>
#include <Geom_Curve.hxx>
#include <Handle_Geom_Curve.hxx>
+#include <Handle_AIS_InteractiveObject.hxx>
+#include <TopoDS_Vertex.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Face.hxx>
+
#include <gp_Pnt2d.hxx>
#include <TopTools_IndexedMapOfShape.hxx>
+#include <TopoDS_Shape.hxx>
+#include <BRep_Builder.hxx>
#include <BRepTools.hxx>
+
#include <TopTools_DataMapOfShapeInteger.hxx>
#include <GProp_GProps.hxx>
#include <BRepGProp.hxx>
#include <fenv.h>
#endif
+#include <Standard_ErrorHandler.hxx>
#include <GeomAPI_ProjectPointOnCurve.hxx>
#include <GeomAPI_ProjectPointOnSurf.hxx>
#include <gp_XY.hxx>
TopTools_IndexedMapOfShape FacesWithEnforcedVertices;
std::map< int, std::set< std::vector<double> > > FaceId2EnforcedVertexCoords;
+std::map< std::vector<double>, std::vector<double> > EnfVertex2ProjVertex;
bool HasSizeMapOnFace=false;
bool HasSizeMapOnEdge=false;
myStudy = NULL;
myStudy = aStudyMgr->GetStudyByID(_studyId);
- MESSAGE("myStudy->StudyId() = " << myStudy->StudyId());
+ if (myStudy)
+ MESSAGE("myStudy->StudyId() = " << myStudy->StudyId());
/* Initialize the Python interpreter */
assert(Py_IsInitialized());
FaceId2AttractorCoords.clear();
FacesWithEnforcedVertices.Clear();
FaceId2EnforcedVertexCoords.clear();
-
+ EnfVertex2ProjVertex.clear();
}
//=============================================================================
double my_u_min=1e6,my_v_min=1e6,my_u_max=-1e6,my_v_max=-1e6;
typedef struct {
- gp_XY uv;
- gp_XYZ xyz;
+ gp_XY uv;
+ gp_XYZ xyz;
} projectionPoint;
/////////////////////////////////////////////////////////
projectionPoint getProjectionPoint(const TopoDS_Face& face, const gp_XYZ& point)
Handle(Geom_Surface) surface = BRep_Tool::Surface(face);
GeomAPI_ProjectPointOnSurf projector( point, surface );
if ( !projector.IsDone() || projector.NbPoints()==0 )
- throw "Can't project";
+ throw "getProjectionPoint: Can't project";
Quantity_Parameter u,v;
projector.LowerDistanceParameters(u,v);
/////////////////////////////////////////////////////////
TopoDS_Shape BLSURFPlugin_BLSURF::entryToShape(std::string entry)
{
- MESSAGE("BLSURFPlugin_BLSURF::entryToShape"<<entry );
+ MESSAGE("BLSURFPlugin_BLSURF::entryToShape "<<entry );
GEOM::GEOM_Object_var aGeomObj;
TopoDS_Shape S = TopoDS_Shape();
SALOMEDS::SObject_var aSObj = myStudy->FindObjectID( entry.c_str() );
}
/////////////////////////////////////////////////////////
-void createEnforcedVertexOnFace(TopoDS_Shape GeomShape, BLSURFPlugin_Hypothesis::TEnforcedVertexList enforcedVertexList)
+void createEnforcedVertexOnFace(TopoDS_Shape GeomShape, BLSURFPlugin_Hypothesis::TEnfVertexList enfVertexList)
{
double xe, ye, ze;
std::vector<double> coords;
- BLSURFPlugin_Hypothesis::TEnforcedVertex enforcedVertex;
- BLSURFPlugin_Hypothesis::TEnforcedVertexList::const_iterator evlIt = enforcedVertexList.begin();
+ std::vector<double> s_coords;
+ std::vector<double> enfVertex;
+// BLSURFPlugin_Hypothesis::TEnfVertex enfVertex;
+ BLSURFPlugin_Hypothesis::TEnfVertexList::const_iterator evlIt = enfVertexList.begin();
- for( ; evlIt != enforcedVertexList.end() ; ++evlIt ) {
+ for( ; evlIt != enfVertexList.end() ; ++evlIt ) {
coords.clear();
- enforcedVertex = *evlIt;
- xe = enforcedVertex[0];
- ye = enforcedVertex[1];
- ze = enforcedVertex[2];
+ s_coords.clear();
+ enfVertex = *evlIt;
+ xe = enfVertex[0];
+ ye = enfVertex[1];
+ ze = enfVertex[2];
MESSAGE("Enforced Vertex: " << xe << ", " << ye << ", " << ze);
// Get the (u,v) values of the enforced vertex on the face
projectionPoint myPoint = getProjectionPoint(TopoDS::Face(GeomShape),gp_XYZ(xe,ye,ze));
Standard_Real y0 = xyzPoint.Y();
Standard_Real z0 = xyzPoint.Z();
MESSAGE("Projected Vertex: " << x0 << ", " << y0 << ", " << z0);
+ MESSAGE("Parametric coordinates: " << u0 << ", " << v0 );
coords.push_back(u0);
coords.push_back(v0);
coords.push_back(x0);
coords.push_back(y0);
coords.push_back(z0);
-
+
+ s_coords.push_back(x0);
+ s_coords.push_back(y0);
+ s_coords.push_back(z0);
+
+ // Save pair projected vertex / enf vertex
+ MESSAGE("Storing pair projected vertex / enf vertex:");
+ MESSAGE("("<< x0 << ", " << y0 << ", " << z0 <<") / (" << xe << ", " << ye << ", " << ze<<")");
+ EnfVertex2ProjVertex[s_coords] = enfVertex;
+
int key = 0;
if (! FacesWithEnforcedVertices.Contains(TopoDS::Face(GeomShape))) {
key = FacesWithEnforcedVertices.Add(TopoDS::Face(GeomShape));
/////////////////////////////////////////////////////////
-void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp, blsurf_session_t *bls)
+void BLSURFPlugin_BLSURF::SetParameters(const BLSURFPlugin_Hypothesis* hyp,
+ blsurf_session_t * bls,
+ SMESH_Mesh& mesh)
{
int _topology = BLSURFPlugin_Hypothesis::GetDefaultTopology();
int _physicalMesh = BLSURFPlugin_Hypothesis::GetDefaultPhysicalMesh();
}
} else {
+ //0020968: EDF1545 SMESH: Problem in the creation of a mesh group on geometry
+ // GetDefaultPhySize() sometimes leads to computation failure
+ _phySize = mesh.GetShapeDiagonalSize() / _gen->GetBoundaryBoxSegmentation();
MESSAGE("BLSURFPlugin_BLSURF::SetParameters using defaults");
}
_smp_phy_size = _phySize;
}
else {
key = VerticesWithSizeMap.FindIndex(TopoDS::Vertex(it.Value()));
-// MESSAGE("Vertex with key " << key << " already in map");
+ MESSAGE("Group of vertices with key " << key << " already in map");
}
+ MESSAGE("Group of vertices with key " << key << " has a size map: " << smIt->second);
VertexId2SizeMap[key] = smIt->second;
}
}
}
else {
key = VerticesWithSizeMap.FindIndex(TopoDS::Vertex(GeomShape));
-// MESSAGE("Vertex with key " << key << " already in map");
+ MESSAGE("Vertex with key " << key << " already in map");
}
+ MESSAGE("Vertex with key " << key << " has a size map: " << smIt->second);
VertexId2SizeMap[key] = smIt->second;
}
}
}
}
}
-
+
if (GeomType == TopAbs_FACE){
HasSizeMapOnFace = true;
createAttractorOnFace(GeomShape, atIt->second);
// Enforced Vertices
//
MESSAGE("Setting Enforced Vertices");
- const BLSURFPlugin_Hypothesis::TEnforcedVertexMap enforcedVertexMap = BLSURFPlugin_Hypothesis::GetAllEnforcedVertices(hyp);
- BLSURFPlugin_Hypothesis::TEnforcedVertexMap::const_iterator enfIt = enforcedVertexMap.begin();
- for ( ; enfIt != enforcedVertexMap.end(); ++enfIt ) {
+ const BLSURFPlugin_Hypothesis::TEntryEnfVertexListMap entryEnfVertexListMap = BLSURFPlugin_Hypothesis::GetAllEnforcedVertices(hyp);
+ BLSURFPlugin_Hypothesis::TEntryEnfVertexListMap::const_iterator enfIt = entryEnfVertexListMap.begin();
+ for ( ; enfIt != entryEnfVertexListMap.end(); ++enfIt ) {
if ( !enfIt->second.empty() ) {
GeomShape = entryToShape(enfIt->first);
GeomType = GeomShape.ShapeType();
status_t curv_fun(real t, real *uv, real *dt, real *dtt, void *user_data);
status_t surf_fun(real *uv, real *xyz, real*du, real *dv,
- real *duu, real *duv, real *dvv, void *user_data);
-status_t message_callback(message_t *msg, void *user_data);
+ real *duu, real *duv, real *dvv, void *user_data);
+status_t message_cb(message_t *msg, void *user_data);
+status_t interrupt_cb(integer *interrupt_status, void *user_data);
//=============================================================================
/*!
MESSAGE("BLSURFPlugin_BLSURF::Compute");
- if (aShape.ShapeType() == TopAbs_COMPOUND) {
- MESSAGE(" the shape is a COMPOUND");
- }
- else {
- MESSAGE(" the shape is UNKNOWN");
- };
+// if (aShape.ShapeType() == TopAbs_COMPOUND) {
+// MESSAGE(" the shape is a COMPOUND");
+// }
+// else {
+// MESSAGE(" the shape is UNKNOWN");
+// };
+ // Fix problem with locales
+ Kernel_Utils::Localizer loc;
+
+ /* create a distene context (generic object) */
+ status_t status = STATUS_ERROR;
+
context_t *ctx = context_new();
- context_set_message_callback(ctx, message_callback, &_comment);
+
+ /* Set the message callback in the working context */
+ context_set_message_callback(ctx, message_cb, &_comment);
+ context_set_interrupt_callback(ctx, interrupt_cb, NULL);
+ /* create the CAD object we will work on. It is associated to the context ctx. */
cad_t *c = cad_new(ctx);
blsurf_session_t *bls = blsurf_session_new(ctx);
VertexId2SizeMap.clear();
MESSAGE("BEGIN SetParameters");
- SetParameters(_hypothesis, bls);
+ SetParameters(_hypothesis, bls, aMesh);
MESSAGE("END SetParameters");
- TopTools_IndexedMapOfShape fmap;
- TopTools_IndexedMapOfShape emap;
- TopTools_IndexedMapOfShape pmap;
+ /* Now fill the CAD object with data from your CAD
+ * environement. This is the most complex part of a successfull
+ * integration.
+ */
+
+ // needed to prevent the opencascade memory managmement from freeing things
vector<Handle(Geom2d_Curve)> curves;
vector<Handle(Geom_Surface)> surfaces;
+ surfaces.resize(0);
+ curves.resize(0);
+
+ TopTools_IndexedMapOfShape fmap;
+ TopTools_IndexedMapOfShape emap;
+ TopTools_IndexedMapOfShape pmap;
+
fmap.Clear();
FaceId2PythonSmp.clear();
emap.Clear();
EdgeId2PythonSmp.clear();
pmap.Clear();
VertexId2PythonSmp.clear();
- surfaces.resize(0);
- curves.resize(0);
assert(Py_IsInitialized());
PyGILState_STATE gstate;
int iface = 0;
string bad_end = "return";
int faceKey = -1;
+ int ienf = 0;
for (TopExp_Explorer face_iter(aShape,TopAbs_FACE);face_iter.More();face_iter.Next()) {
TopoDS_Face f=TopoDS::Face(face_iter.Current());
+
+ // make INTERNAL face oriented FORWARD (issue 0020993)
+ if (f.Orientation() != TopAbs_FORWARD && f.Orientation() != TopAbs_REVERSED )
+ f.Orientation(TopAbs_FORWARD);
if (fmap.FindIndex(f) > 0)
continue;
fmap.Add(f);
iface++;
surfaces.push_back(BRep_Tool::Surface(f));
-
+
+ /* create an object representing the face for blsurf */
+ /* where face_id is an integer identifying the face.
+ * surf_function is the function that defines the surface
+ * (For this face, it will be called by blsurf with your_face_object_ptr
+ * as last parameter.
+ */
cad_face_t *fce = cad_face_new(c, iface, surf_fun, surfaces.back());
+
+ /* by default a face has no tag (color). The following call sets it to the same value as the face_id : */
cad_face_set_tag(fce, iface);
+
+ /* Set face orientation (optional if you want a well oriented output mesh)*/
if(f.Orientation() != TopAbs_FORWARD){
cad_face_set_orientation(fce, CAD_ORIENTATION_REVERSED);
} else {
}
if (HasSizeMapOnFace){
- std::cout << "A size map is defined on a face" << std::endl;
+// std::cout << "A size map is defined on a face" << std::endl;
// Classic size map
faceKey = FacesWithSizeMap.FindIndex(f);
faceKey = FacesWithEnforcedVertices.FindIndex(f);
std::map<int,std::set<std::vector<double> > >::const_iterator evmIt = FaceId2EnforcedVertexCoords.find(faceKey);
if (evmIt != FaceId2EnforcedVertexCoords.end()) {
- std::cout << "Some enforced vertices are defined" << std::endl;
- int ienf = 0;
+ MESSAGE("Some enforced vertices are defined");
+// int ienf = 0;
std::set<std::vector<double> > evl;
// std::vector<double> ev;
MESSAGE("Face indice: " << iface);
MESSAGE("Adding enforced vertices");
evl = evmIt->second;
- MESSAGE("Number of vertices to add: "<< evl.size())
+ MESSAGE("Number of vertices to add: "<< evl.size());
std::set< std::vector<double> >::const_iterator evlIt = evl.begin();
for (; evlIt != evl.end(); ++evlIt) {
-// ev = *evlIt;
-// for (int i=0; i<evl.size() ; i++) {
-// ev = evl[i];
-
-// double xyzCoords[3] = {ev[2], ev[3], ev[4]};
- double xyzCoords[3] = {evlIt->at(0), evlIt->at(3), evlIt->at(4)};
+ std::vector<double> xyzCoords;
+ xyzCoords.push_back(evlIt->at(2));
+ xyzCoords.push_back(evlIt->at(3));
+ xyzCoords.push_back(evlIt->at(4));
+// double xyzCoords[3] = {evlIt->at(2), evlIt->at(3), evlIt->at(4)};
MESSAGE("Check position of vertex =(" << xyzCoords[0] << "," << xyzCoords[1] << "," << xyzCoords[2] << ")");
gp_Pnt P(xyzCoords[0],xyzCoords[1],xyzCoords[2]);
BRepClass_FaceClassifier scl(f,P,1e-7);
BRepClass_FaceClassifierPerform(&scl,f,P,1e-7);
TopAbs_State result = scl.State();
MESSAGE("Position of point on face: "<<result);
- if ( result == TopAbs_OUT )
+ if ( result == TopAbs_OUT ) {
MESSAGE("Point is out of face: node is not created");
- if ( result == TopAbs_UNKNOWN )
- MESSAGE("Point position on face is unknown: node is not created");
- if ( result == TopAbs_ON )
- MESSAGE("Point is on border of face: node is not created");
+ if (EnfVertex2ProjVertex.find(xyzCoords) != EnfVertex2ProjVertex.end())
+ EnfVertex2ProjVertex.erase(xyzCoords);
+ }
+ if ( result == TopAbs_UNKNOWN ) {
+ MESSAGE("Point position on face is unknown: node is not created");
+ if (EnfVertex2ProjVertex.find(xyzCoords) != EnfVertex2ProjVertex.end())
+ EnfVertex2ProjVertex.erase(xyzCoords);
+ }
+ if ( result == TopAbs_ON ) {
+ MESSAGE("Point is on border of face: node is not created");
+ if (EnfVertex2ProjVertex.find(xyzCoords) != EnfVertex2ProjVertex.end())
+ EnfVertex2ProjVertex.erase(xyzCoords);
+ }
if ( result == TopAbs_IN )
{
// Point is inside face and not on border
}
FaceId2EnforcedVertexCoords.erase(faceKey);
}
- else
- std::cout << "No enforced vertex defined" << std::endl;
+// else
+// std::cout << "No enforced vertex defined" << std::endl;
}
/****************************************************************************************
EDGES
+ now create the edges associated to this face
*****************************************************************************************/
int edgeKey = -1;
for (TopExp_Explorer edge_iter(f,TopAbs_EDGE);edge_iter.More();edge_iter.Next()) {
if (HasSizeMapOnEdge){
edgeKey = EdgesWithSizeMap.FindIndex(e);
if (EdgeId2SizeMap.find(edgeKey)!=EdgeId2SizeMap.end()) {
- theSizeMapStr = EdgeId2SizeMap[faceKey];
+ MESSAGE("Sizemap defined on edge with index " << edgeKey);
+ theSizeMapStr = EdgeId2SizeMap[edgeKey];
if (theSizeMapStr.find(bad_end) == (theSizeMapStr.size()-bad_end.size()-1))
continue;
// Expr To Python function, verification is performed at validation in GUI
EdgeId2SizeMap.erase(edgeKey);
}
}
+
+ /* attach the edge to the current blsurf face */
cad_edge_t *edg = cad_edge_new(fce, ic, tmin, tmax, curv_fun, curves.back());
+
+ /* by default an edge has no tag (color). The following call sets it to the same value as the edge_id : */
cad_edge_set_tag(edg, ic);
+
+ /* by default, an edge does not necessalry appear in the resulting mesh,
+ unless the following property is set :
+ */
cad_edge_set_property(edg, EDGE_PROPERTY_SOFT_REQUIRED);
+
+ /* by default an edge is a boundary edge */
if (e.Orientation() == TopAbs_INTERNAL)
cad_edge_set_property(edg, EDGE_PROPERTY_INTERNAL);
if(*ip <= 0)
*ip = pmap.Add(v);
- vertexKey = VerticesWithSizeMap.FindIndex(v);
+ //vertexKey = VerticesWithSizeMap.FindIndex(v);
if (HasSizeMapOnVertex){
vertexKey = VerticesWithSizeMap.FindIndex(v);
if (VertexId2SizeMap.find(vertexKey)!=VertexId2SizeMap.end()){
- theSizeMapStr = VertexId2SizeMap[faceKey];
+ theSizeMapStr = VertexId2SizeMap[vertexKey];
+ //MESSAGE("VertexId2SizeMap[faceKey]: " << VertexId2SizeMap[vertexKey]);
if (theSizeMapStr.find(bad_end) == (theSizeMapStr.size()-bad_end.size()-1))
continue;
// Expr To Python function, verification is performed at validation in GUI
PyObject * func = NULL;
func = PyObject_GetAttrString(main_mod, "f");
VertexId2PythonSmp[*ip]=func;
-// VertexId2SizeMap.erase(vertexKey); // do not erase if using a vector
+ VertexId2SizeMap.erase(vertexKey); // do not erase if using a vector
}
}
}
// should not happen
MESSAGE("An edge does not have 2 extremities.");
} else {
- if (d1 < d2)
+ if (d1 < d2) {
+ // This defines the curves extremity connectivity
cad_edge_set_extremities(edg, ip1, ip2);
- else
+ /* set the tag (color) to the same value as the extremity id : */
+ cad_edge_set_extremities_tag(edg, ip1, ip2);
+ }
+ else {
cad_edge_set_extremities(edg, ip2, ip1);
+ cad_edge_set_extremities_tag(edg, ip2, ip1);
+ }
}
} // for edge
} //for face
int oldFEFlags = fedisableexcept( FE_ALL_EXCEPT );
#endif
- status_t status = STATUS_ERROR;
-
try {
OCC_CATCH_SIGNALS;
_comment = "Exception in blsurf_compute_mesh()";
}
if ( status != STATUS_OK) {
+ // There was an error while meshing
blsurf_session_delete(bls);
cad_delete(c);
context_delete(ctx);
std::cout << "End of Surface Mesh generation" << std::endl;
std::cout << std::endl;
- mesh_t *msh;
+ mesh_t *msh = NULL;
blsurf_data_get_mesh(bls, &msh);
if(!msh){
blsurf_session_delete(bls);
//return false;
}
+ /* retrieve mesh data (see distene/mesh.h) */
integer nv, ne, nt, nq, vtx[4], tag;
real xyz[3];
SMDS_MeshNode** nodes = new SMDS_MeshNode*[nv+1];
bool* tags = new bool[nv+1];
+ /* enumerated vertices */
for(int iv=1;iv<=nv;iv++) {
mesh_get_vertex_coordinates(msh, iv, xyz);
mesh_get_vertex_tag(msh, iv, &tag);
+ // Issue 0020656. Use vertex coordinates
+ if ( tag > 0 && tag <= pmap.Extent() ) {
+ TopoDS_Vertex v = TopoDS::Vertex(pmap(tag));
+ double tol = BRep_Tool::Tolerance( v );
+ gp_Pnt p = BRep_Tool::Pnt( v );
+ if ( p.IsEqual( gp_Pnt( xyz[0], xyz[1], xyz[2]), 2*tol))
+ xyz[0] = p.X(), xyz[1] = p.Y(), xyz[2] = p.Z();
+ else
+ tag = 0; // enforced or attracted vertex
+ }
nodes[iv] = meshDS->AddNode(xyz[0], xyz[1], xyz[2]);
+
+ /* TODO GROUPS
+ // Create group of enforced vertices if requested
+ if(_hypothesis) {
+ std::vector<double> projVertex;
+ projVertex.push_back((double)xyz[0]);
+ projVertex.push_back((double)xyz[1]);
+ projVertex.push_back((double)xyz[2]);
+ std::map< std::vector<double>, std::vector<double> >::const_iterator projIt = EnfVertex2ProjVertex.find(projVertex);
+ if (projIt != EnfVertex2ProjVertex.end()) {
+ double x = (projIt->second)[0];
+ double y = (projIt->second)[1];
+ double z = (projIt->second)[2];
+ BLSURFPlugin_Hypothesis::TEnfVertex enfVertex;
+ enfVertex.push_back(x);
+ enfVertex.push_back(y);
+ enfVertex.push_back(z);
+
+ BLSURFPlugin_Hypothesis::TEnfVertexGroupNameMap groupNameMap = _hypothesis->_GetEnforcedVertexGroupNameMap();
+ BLSURFPlugin_Hypothesis::TEnfVertexGroupNameMap::const_iterator groupNameMapIt = groupNameMap.find(enfVertex);
+ if (groupNameMapIt != groupNameMap.end()) {
+ MESSAGE("Found enforced vertex @ " << xyz[0] << ", " << xyz[1] << ", " << xyz[2])
+ BLSURFPlugin_Hypothesis::TEnfGroupName groupName = groupNameMapIt->second;
+ if (groupName != "") {
+ bool groupDone = false;
+ const set<SMESHDS_GroupBase*>& allGroups = meshDS->GetGroups();
+ set<SMESHDS_GroupBase*>::const_iterator grIt;
+ MESSAGE("Parsing the groups of the mesh");
+ for ( grIt = allGroups.begin(); grIt != allGroups.end(); ++grIt ) {
+ SMESHDS_Group* group = dynamic_cast<SMESHDS_Group*>( *grIt );
+ MESSAGE("Group: " << group->GetStoreName());
+ if ( group && group->SMDSGroup().GetType()==SMDSAbs_Node
+ && groupName.compare(group->GetStoreName())==0) {
+ group->SMDSGroup().Add(nodes[iv]);
+// int id = // recuperer l'id SMESH du noeud
+// _hypothesis->AddEnfVertexIDs(groupName,id)
+ groupDone = true;
+ MESSAGE("Successfully added enforced vertex to existing group " << groupName);
+ break;
+ }
+ }
+ if (!groupDone)
+ {
+ int groupId;
+ SMESH_Group* aGroup = aMesh.AddGroup(SMDSAbs_Node, groupName.c_str(), groupId);
+ if ( aGroup ) {
+ SMESHDS_Group* aGroupDS = dynamic_cast<SMESHDS_Group*>( aGroup->GetGroupDS() );
+ if ( aGroupDS ) {
+ aGroupDS->SetStoreName( groupName.c_str() );
+ aGroupDS->SMDSGroup().Add(nodes[iv]);
+ MESSAGE("Successfully created enforced vertex group " << groupName);
+ groupDone = true;
+ }
+ }
+ }
+ if (!groupDone)
+ throw SALOME_Exception(LOCALIZED("A enforced vertex node was not added to a group"));
+ }
+ else
+ MESSAGE("Group name is empty: '"<<groupName<<"' => group is not created");
+ }
+ else
+ MESSAGE("No group name for projected vertex ("<<x<<","<<y<<","<<z<<")")
+ }
+// else
+// MESSAGE("No group name for vertex ("<<xyz[0]<<","<<xyz[1]<<","<<xyz[2]<<")")
+ }
+ */
+
// internal point are tagged to zero
- if(tag){
+ if(tag > 0 && tag <= pmap.Extent() ){
meshDS->SetNodeOnVertex(nodes[iv], TopoDS::Vertex(pmap(tag)));
tags[iv] = false;
} else {
}
}
+ /* enumerate edges */
for(int it=1;it<=ne;it++) {
mesh_get_edge_vertices(msh, it, vtx);
SMDS_MeshEdge* edg = meshDS->AddEdge(nodes[vtx[0]], nodes[vtx[1]]);
}
+ /* enumerate triangles */
for(int it=1;it<=nt;it++) {
mesh_get_triangle_vertices(msh, it, vtx);
SMDS_MeshFace* tri = meshDS->AddFace(nodes[vtx[0]], nodes[vtx[1]], nodes[vtx[2]]);
};
}
+ /* enumerate quadrangles */
for(int it=1;it<=nq;it++) {
mesh_get_quadrangle_vertices(msh, it, vtx);
SMDS_MeshFace* quad = meshDS->AddFace(nodes[vtx[0]], nodes[vtx[1]], nodes[vtx[2]], nodes[vtx[3]]);
feenableexcept( oldFEFlags );
feclearexcept( FE_ALL_EXCEPT );
#endif
-
+
+ /*
std::cout << "FacesWithSizeMap" << std::endl;
FacesWithSizeMap.Statistics(std::cout);
std::cout << "EdgesWithSizeMap" << std::endl;
VerticesWithSizeMap.Statistics(std::cout);
std::cout << "FacesWithEnforcedVertices" << std::endl;
FacesWithEnforcedVertices.Statistics(std::cout);
+ */
return true;
}
Standard_Real p0 = 0.0;
Standard_Real p1 = 1.0;
- Handle(Geom_Curve) curve = BRep_Tool::Curve(edge, p0, p1);
-
- GeomAPI_ProjectPointOnCurve proj(pnt, curve);
+ TopLoc_Location loc;
+ Handle(Geom_Curve) curve = BRep_Tool::Curve(edge, loc, p0, p1);
- double pa = (double)proj.Parameter(1);
+ if ( !loc.IsIdentity() ) pnt.Transform( loc.Transformation().Inverted() );
+ GeomAPI_ProjectPointOnCurve proj(pnt, curve, p0, p1);
- GProp_GProps LProps;
- BRepGProp::LinearProperties(ed, LProps);
- double lg = (double)LProps.Mass();
+ double pa = 0.;
+ if ( proj.NbPoints() > 0 )
+ {
+ pa = (double)proj.LowerDistanceParameter();
+ // Issue 0020656. Move node if it is too far from edge
+ gp_Pnt curve_pnt = curve->Value( pa );
+ double dist2 = pnt.SquareDistance( curve_pnt );
+ double tol = BRep_Tool::Tolerance( edge );
+ if ( 1e-12 < dist2 && dist2 <= 2*tol*tol ) // large enough and within tolerance
+ {
+ curve_pnt.Transform( loc );
+ meshDS->MoveNode( node, curve_pnt.X(), curve_pnt.Y(), curve_pnt.Z() );
+ }
+ }
+// GProp_GProps LProps;
+// BRepGProp::LinearProperties(ed, LProps);
+// double lg = (double)LProps.Mass();
meshDS->SetNodeOnEdge(node, edge, pa);
}
return hyp.LoadFrom( load );
}
+/* Curve definition function See cad_curv_t in file distene/cad.h for
+ * more information.
+ * NOTE : if when your CAD systems evaluates second
+ * order derivatives it also computes first order derivatives and
+ * function evaluation, you can optimize this example by making only
+ * one CAD call and filling the necessary uv, dt, dtt arrays.
+ */
status_t curv_fun(real t, real *uv, real *dt, real *dtt, void *user_data)
{
+ /* t is given. It contains the t (time) 1D parametric coordintaes
+ of the point PreCAD/BLSurf is querying on the curve */
+
+ /* user_data identifies the edge PreCAD/BLSurf is querying
+ * (see cad_edge_new later in this example) */
const Geom2d_Curve*pargeo = (const Geom2d_Curve*) user_data;
if (uv){
+ /* BLSurf is querying the function evaluation */
gp_Pnt2d P;
P=pargeo->Value(t);
uv[0]=P.X(); uv[1]=P.Y();
}
if(dt) {
+ /* query for the first order derivatives */
gp_Vec2d V1;
V1=pargeo->DN(t,1);
dt[0]=V1.X(); dt[1]=V1.Y();
}
if(dtt){
+ /* query for the second order derivatives */
gp_Vec2d V2;
V2=pargeo->DN(t,2);
dtt[0]=V2.X(); dtt[1]=V2.Y();
}
- return 0;
+ return STATUS_OK;
}
+/* Surface definition function.
+ * See cad_surf_t in file distene/cad.h for more information.
+ * NOTE : if when your CAD systems evaluates second order derivatives it also
+ * computes first order derivatives and function evaluation, you can optimize
+ * this example by making only one CAD call and filling the necessary xyz, du, dv, etc..
+ * arrays.
+ */
status_t surf_fun(real *uv, real *xyz, real*du, real *dv,
- real *duu, real *duv, real *dvv, void *user_data)
+ real *duu, real *duv, real *dvv, void *user_data)
{
+ /* uv[2] is given. It contains the u,v coordinates of the point
+ * PreCAD/BLSurf is querying on the surface */
+
+ /* user_data identifies the face PreCAD/BLSurf is querying (see
+ * cad_face_new later in this example)*/
const Geom_Surface* geometry = (const Geom_Surface*) user_data;
if(xyz){
}
if(duu && duv && dvv){
+
gp_Pnt P;
gp_Vec D1U,D1V;
gp_Vec D2U,D2V,D2UV;
dvv[0]=D2V.X(); dvv[1]=D2V.Y(); dvv[2]=D2V.Z();
}
- return 0;
+ return STATUS_OK;
}
return STATUS_OK;
}
-status_t message_callback(message_t *msg, void *user_data)
+/*
+ * The following function will be called for PreCAD/BLSurf message
+ * printing. See context_set_message_callback (later in this
+ * template) for how to set user_data.
+ */
+status_t message_cb(message_t *msg, void *user_data)
{
integer errnumber = 0;
char *desc;
return STATUS_OK;
}
+/* This is the interrupt callback. PreCAD/BLSurf will call this
+ * function regularily. See the file distene/interrupt.h
+ */
+status_t interrupt_cb(integer *interrupt_status, void *user_data)
+{
+ integer you_want_to_continue = 1;
+
+ if(you_want_to_continue)
+ *interrupt_status = INTERRUPT_CONTINUE;
+ else /* you want to stop BLSurf */
+ *interrupt_status = INTERRUPT_STOP;
+
+ return STATUS_OK;
+}
//=============================================================================
/*!
*
*/
//=============================================================================
-bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh& aMesh,
- const TopoDS_Shape& aShape,
- MapShapeNbElems& aResMap)
+bool BLSURFPlugin_BLSURF::Evaluate(SMESH_Mesh& aMesh,
+ const TopoDS_Shape& aShape,
+ MapShapeNbElems& aResMap)
{
int _physicalMesh = BLSURFPlugin_Hypothesis::GetDefaultPhysicalMesh();
double _phySize = BLSURFPlugin_Hypothesis::GetDefaultPhySize();
//int _geometricMesh = BLSURFPlugin_Hypothesis::GetDefaultGeometricMesh();
//double _angleMeshS = BLSURFPlugin_Hypothesis::GetDefaultAngleMeshS();
double _angleMeshC = BLSURFPlugin_Hypothesis::GetDefaultAngleMeshC();
+ bool _quadAllowed = BLSURFPlugin_Hypothesis::GetDefaultQuadAllowed();
if(_hypothesis) {
_physicalMesh = (int) _hypothesis->GetPhysicalMesh();
_phySize = _hypothesis->GetPhySize();
//_geometricMesh = (int) hyp->GetGeometricMesh();
//_angleMeshS = hyp->GetAngleMeshS();
_angleMeshC = _hypothesis->GetAngleMeshC();
+ _quadAllowed = _hypothesis->GetQuadAllowed();
}
bool IsQuadratic = false;
C->D0(f+dp,P2);
gp_Vec V1(P1,P2);
for(int j=2; j<=200; j++) {
- C->D0(f+dp*j,P3);
- gp_Vec V2(P2,P3);
- fullAng += fabs(V1.Angle(V2));
- V1 = V2;
- P2 = P3;
+ C->D0(f+dp*j,P3);
+ gp_Vec V2(P2,P3);
+ fullAng += fabs(V1.Angle(V2));
+ V1 = V2;
+ P2 = P3;
}
nb1d = (int)( fullAng/_angleMeshC + 1 );
}
BRepGProp::SurfaceProperties(F,G);
double anArea = G.Mass();
int nb1d = 0;
+ std::vector<int> nb1dVec;
for (TopExp_Explorer exp1(F,TopAbs_EDGE); exp1.More(); exp1.Next()) {
- nb1d += EdgesMap.Find(exp1.Current());
+ int nbSeg = EdgesMap.Find(exp1.Current());
+ nb1d += nbSeg;
+ nb1dVec.push_back( nbSeg );
}
- int nbFaces = (int) ( anArea/( ELen*ELen*sqrt(3.) / 4 ) );
- int nbNodes = (int) ( ( nbFaces*3 - (nb1d-1)*2 ) / 6 + 1 );
- std::vector<int> aVec(SMDSEntity_Last);
- for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
+ int nbQuad = 0;
+ int nbTria = (int) ( anArea/( ELen*ELen*sqrt(3.) / 4 ) );
+ int nbNodes = (int) ( ( nbTria*3 - (nb1d-1)*2 ) / 6 + 1 );
+ if ( _quadAllowed )
+ {
+ if ( nb1dVec.size() == 4 ) // quadrangle geom face
+ {
+ int n1 = nb1dVec[0], n2 = nb1dVec[ nb1dVec[1] == nb1dVec[0] ? 2 : 1 ];
+ nbQuad = n1 * n2;
+ nbNodes = (n1 + 1) * (n2 + 1);
+ nbTria = 0;
+ }
+ else
+ {
+ nbTria = nbQuad = nbTria / 3 + 1;
+ }
+ }
+ std::vector<int> aVec(SMDSEntity_Last,0);
if( IsQuadratic ) {
- int nb1d_in = (nbFaces*3 - nb1d) / 2;
+ int nb1d_in = (nbTria*3 - nb1d) / 2;
aVec[SMDSEntity_Node] = nbNodes + nb1d_in;
- aVec[SMDSEntity_Quad_Triangle] = nbFaces;
+ aVec[SMDSEntity_Quad_Triangle] = nbTria;
+ aVec[SMDSEntity_Quad_Quadrangle] = nbQuad;
}
else {
aVec[SMDSEntity_Node] = nbNodes;
- aVec[SMDSEntity_Triangle] = nbFaces;
+ aVec[SMDSEntity_Triangle] = nbTria;
+ aVec[SMDSEntity_Quadrangle] = nbQuad;
}
aResMap.insert(std::make_pair(sm,aVec));
}
BRepGProp::VolumeProperties(aShape,G);
double aVolume = G.Mass();
double tetrVol = 0.1179*ELen*ELen*ELen;
- int nbVols = (int)aVolume/tetrVol;
- int nb1d_in = (int) ( ( nbVols*6 - fullNbSeg ) / 6 );
+ int nbVols = int(aVolume/tetrVol);
+ int nb1d_in = int(( nbVols*6 - fullNbSeg ) / 6 );
std::vector<int> aVec(SMDSEntity_Last);
for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
if( IsQuadratic ) {
