From: rnv Date: Tue, 30 May 2017 12:58:50 +0000 (+0300) Subject: Merge rnv/win_swig_generation branch. X-Git-Tag: Before_multi_study_removal_06072017~1 X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=commitdiff_plain;h=d92118a7a570fd784e024a7a893a67a4fc8f112c;hp=82e470f3c19548107c7f19ad1641dcfc368fc313 Merge rnv/win_swig_generation branch. --- diff --git a/CMakeLists.txt b/CMakeLists.txt index ef42ec11a..b3e766d8d 100755 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -35,7 +35,7 @@ SET(${PROJECT_NAME_UC}_MINOR_VERSION 3) SET(${PROJECT_NAME_UC}_PATCH_VERSION 0) SET(${PROJECT_NAME_UC}_VERSION ${${PROJECT_NAME_UC}_MAJOR_VERSION}.${${PROJECT_NAME_UC}_MINOR_VERSION}.${${PROJECT_NAME_UC}_PATCH_VERSION}) -SET(${PROJECT_NAME_UC}_VERSION_DEV 0) +SET(${PROJECT_NAME_UC}_VERSION_DEV 1) # Common CMake macros # =================== diff --git a/doc/salome/gui/SMESH/images/prism_mesh.png b/doc/salome/gui/SMESH/images/prism_mesh.png new file mode 100644 index 000000000..95a31214c Binary files /dev/null and b/doc/salome/gui/SMESH/images/prism_mesh.png differ diff --git a/doc/salome/gui/SMESH/input/constructing_meshes.doc b/doc/salome/gui/SMESH/input/constructing_meshes.doc index 6794b398b..2481807fd 100644 --- a/doc/salome/gui/SMESH/input/constructing_meshes.doc +++ b/doc/salome/gui/SMESH/input/constructing_meshes.doc @@ -19,7 +19,7 @@ and hypotheses. Mesh generation on the geometry is performed in the bottom-up flow: nodes on vertices are created first, then edges are divided into -segments using nodes on vertices; the node of segments are then +segments using nodes on vertices; the nodes of segments are then used to mesh faces; then the nodes of faces are used to mesh solids. This automatically assures the conformity of the mesh. diff --git a/doc/salome/gui/SMESH/input/prism_3d_algo.doc b/doc/salome/gui/SMESH/input/prism_3d_algo.doc index 59aee72fa..6d4be7be8 100644 --- a/doc/salome/gui/SMESH/input/prism_3d_algo.doc +++ b/doc/salome/gui/SMESH/input/prism_3d_algo.doc @@ -6,6 +6,8 @@ Extrusion 3D algorithm can be used for meshing prisms, i.e. 3D shapes defined by two opposing faces having the same number of vertices and edges. These two faces should be connected by quadrangle "side" faces. +\image html prism_mesh.png "Clipping view of a mesh of a prism with non-planar base and top faces" + The prism is allowed to have sides composed of several faces. (A prism side is a row of faces (or one face) connecting the corresponding edges of the top and base faces). However, a prism diff --git a/doc/salome/gui/SMESH/input/smeshpy_interface.doc b/doc/salome/gui/SMESH/input/smeshpy_interface.doc index 1dc47c407..524a4ceb0 100644 --- a/doc/salome/gui/SMESH/input/smeshpy_interface.doc +++ b/doc/salome/gui/SMESH/input/smeshpy_interface.doc @@ -2,10 +2,10 @@ \page smeshpy_interface_page Python interface -Python API for SALOME %Mesh module defines several classes that can +Python API of SALOME %Mesh module defines several classes that can be used for easy mesh creation and edition. -Documentation for SALOME %Mesh module Python API is available in two forms: +Documentation of SALOME %Mesh module Python API is available in two forms: - Structured documentation, where all methods and classes are grouped by their functionality. - Linear documentation grouped only by classes, declared diff --git a/idl/SMESH_Mesh.idl b/idl/SMESH_Mesh.idl index 805b21df6..cef3831ee 100644 --- a/idl/SMESH_Mesh.idl +++ b/idl/SMESH_Mesh.idl @@ -975,6 +975,11 @@ module SMESH */ long FindElementByNodes(in long_array nodes); + /*! + * Return elements including all given nodes. + */ + long_array GetElementsByNodes(in long_array nodes, in ElementType elem_type); + /*! * Returns true if given element is polygon */ diff --git a/src/OBJECT/SMESH_ScalarBarActor.cxx b/src/OBJECT/SMESH_ScalarBarActor.cxx index a45696d7e..8096b1609 100644 --- a/src/OBJECT/SMESH_ScalarBarActor.cxx +++ b/src/OBJECT/SMESH_ScalarBarActor.cxx @@ -28,16 +28,16 @@ #include #include +#include #include #include #include +#include #include #include #include #include #include -#include -#include #define SHRINK_COEF 0.08; @@ -51,13 +51,14 @@ vtkCxxSetObjectMacro(SMESH_ScalarBarActor,TitleTextProperty,vtkTextProperty); // Instantiate object with 64 maximum colors; 5 labels; %%-#6.3g label // format, no title, and vertical orientation. The initial scalar bar // size is (0.05 x 0.8) of the viewport size. -SMESH_ScalarBarActor::SMESH_ScalarBarActor() { +SMESH_ScalarBarActor::SMESH_ScalarBarActor() +{ this->LookupTable = NULL; this->Position2Coordinate->SetValue(0.17, 0.8); - + this->PositionCoordinate->SetCoordinateSystemToNormalizedViewport(); this->PositionCoordinate->SetValue(0.82,0.1); - + this->MaximumNumberOfColors = 64; this->NumberOfLabels = 5; this->NumberOfLabelsBuilt = 0; @@ -74,7 +75,7 @@ SMESH_ScalarBarActor::SMESH_ScalarBarActor() { this->TitleTextProperty = vtkTextProperty::New(); this->TitleTextProperty->ShallowCopy(this->LabelTextProperty); - this->LabelFormat = new char[8]; + this->LabelFormat = new char[8]; sprintf(this->LabelFormat,"%s","%-#6.3g"); this->TitleMapper = vtkTextMapper::New(); @@ -82,7 +83,7 @@ SMESH_ScalarBarActor::SMESH_ScalarBarActor() { this->TitleActor->SetMapper(this->TitleMapper); this->TitleActor->GetPositionCoordinate()-> SetReferenceCoordinate(this->PositionCoordinate); - + this->TextMappers = NULL; this->TextActors = NULL; @@ -104,7 +105,7 @@ SMESH_ScalarBarActor::SMESH_ScalarBarActor() { myDistribution = vtkPolyData::New(); myDistributionMapper = vtkPolyDataMapper2D::New(); myDistributionMapper->SetInputData(this->myDistribution); - + myDistributionActor = vtkActor2D::New(); myDistributionActor->SetMapper(this->myDistributionMapper); myDistributionActor->GetPositionCoordinate()-> @@ -129,12 +130,12 @@ void SMESH_ScalarBarActor::ReleaseGraphicsResources(vtkWindow *win) { this->TitleActor->ReleaseGraphicsResources(win); if (this->TextMappers != NULL ) - { + { for (int i=0; i < this->NumberOfLabelsBuilt; i++) - { + { this->TextActors[i]->ReleaseGraphicsResources(win); - } } + } this->ScalarBarActor->ReleaseGraphicsResources(win); // rnv begin // Customization of the vtkScalarBarActor to show distribution histogram. @@ -143,41 +144,42 @@ void SMESH_ScalarBarActor::ReleaseGraphicsResources(vtkWindow *win) /*--------------------------------------------------------------------------*/ -SMESH_ScalarBarActor::~SMESH_ScalarBarActor() { - if (this->LabelFormat) - { +SMESH_ScalarBarActor::~SMESH_ScalarBarActor() +{ + if (this->LabelFormat) + { delete [] this->LabelFormat; this->LabelFormat = NULL; - } + } this->TitleMapper->Delete(); this->TitleActor->Delete(); if (this->TextMappers != NULL ) - { + { for (int i=0; i < this->NumberOfLabelsBuilt; i++) - { + { this->TextMappers[i]->Delete(); this->TextActors[i]->Delete(); - } + } delete [] this->TextMappers; delete [] this->TextActors; - } + } this->ScalarBar->Delete(); this->ScalarBarMapper->Delete(); this->ScalarBarActor->Delete(); if (this->Title) - { + { delete [] this->Title; this->Title = NULL; - } - + } + this->SetLookupTable(NULL); this->SetLabelTextProperty(NULL); this->SetTitleTextProperty(NULL); - + // rnv begin // Customization of the vtkScalarBarActor to show distribution histogram: myDistribution->Delete(); @@ -191,26 +193,26 @@ int SMESH_ScalarBarActor::RenderOverlay(vtkViewport *viewport) { int renderedSomething = 0; int i; - + // Everything is built, just have to render if (this->Title != NULL) - { + { renderedSomething += this->TitleActor->RenderOverlay(viewport); - } - if (!myTitleOnlyVisibility) { + } + if ( !myTitleOnlyVisibility ) { this->ScalarBarActor->RenderOverlay(viewport); this->myDistributionActor->RenderOverlay(viewport); - if( this->TextActors == NULL) - { - vtkWarningMacro(<<"Need a mapper to render a scalar bar"); - return renderedSomething; - } - - for (i=0; iNumberOfLabels; i++) - { + if ( this->TextActors == NULL ) + { + vtkWarningMacro(<<"Need a mapper to render a scalar bar"); + return renderedSomething; + } + + for ( i=0; iNumberOfLabels; i++ ) + { renderedSomething += this->TextActors[i]->RenderOverlay(viewport); - } - } + } + } renderedSomething = (renderedSomething > 0)?(1):(0); return renderedSomething; @@ -223,75 +225,75 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) int renderedSomething = 0; int i; int size[2]; - + if (!this->LookupTable) - { + { vtkWarningMacro(<<"Need a mapper to render a scalar bar"); return 0; - } + } if (!this->TitleTextProperty) - { + { vtkErrorMacro(<<"Need title text property to render a scalar bar"); return 0; - } + } if (!this->LabelTextProperty) - { + { vtkErrorMacro(<<"Need label text property to render a scalar bar"); return 0; - } + } // Check to see whether we have to rebuild everything int positionsHaveChanged = 0; - if (viewport->GetMTime() > this->BuildTime || - (viewport->GetVTKWindow() && + if (viewport->GetMTime() > this->BuildTime || + (viewport->GetVTKWindow() && viewport->GetVTKWindow()->GetMTime() > this->BuildTime)) - { + { // if the viewport has changed we may - or may not need // to rebuild, it depends on if the projected coords chage int *barOrigin; barOrigin = this->PositionCoordinate->GetComputedViewportValue(viewport); - size[0] = + size[0] = this->Position2Coordinate->GetComputedViewportValue(viewport)[0] - barOrigin[0]; - size[1] = + size[1] = this->Position2Coordinate->GetComputedViewportValue(viewport)[1] - barOrigin[1]; - if (this->LastSize[0] != size[0] || + if (this->LastSize[0] != size[0] || this->LastSize[1] != size[1] || - this->LastOrigin[0] != barOrigin[0] || + this->LastOrigin[0] != barOrigin[0] || this->LastOrigin[1] != barOrigin[1]) - { + { positionsHaveChanged = 1; - } } - + } + // Check to see whether we have to rebuild everything - if (positionsHaveChanged || - this->GetMTime() > this->BuildTime || - this->LookupTable->GetMTime() > this->BuildTime || - this->LabelTextProperty->GetMTime() > this->BuildTime || - this->TitleTextProperty->GetMTime() > this->BuildTime) - { + if ( positionsHaveChanged || + this->GetMTime() > this->BuildTime || + this->LookupTable->GetMTime() > this->BuildTime || + this->LabelTextProperty->GetMTime() > this->BuildTime || + this->TitleTextProperty->GetMTime() > this->BuildTime) + { vtkDebugMacro(<<"Rebuilding subobjects"); // Delete previously constructed objects // - if (this->TextMappers != NULL ) + if ( this->TextMappers != NULL ) + { + for ( i = 0; i < this->NumberOfLabelsBuilt; i++ ) { - for (i=0; i < this->NumberOfLabelsBuilt; i++) - { this->TextMappers[i]->Delete(); this->TextActors[i]->Delete(); - } + } delete [] this->TextMappers; delete [] this->TextActors; - } + } // Build scalar bar object; determine its type // - // is this a vtkLookupTable or a subclass of vtkLookupTable + // is this a vtkLookupTable or a subclass of vtkLookupTable // with its scale set to log // NOTE: it's possible we could to without the 'lut' variable // later in the code, but if the vtkLookupTableSafeDownCast operation @@ -300,13 +302,13 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) vtkLookupTable *LUT = vtkLookupTable::SafeDownCast( this->LookupTable ); int isLogTable = 0; if ( LUT ) - { + { if ( LUT->GetScale() == VTK_SCALE_LOG10 ) - { - isLogTable = 1; - } + { + isLogTable = 1; } - + } + // we hard code how many steps to display vtkScalarsToColors *lut = this->LookupTable; int numColors = this->MaximumNumberOfColors; @@ -332,7 +334,8 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) if(!distrVisibility) vtkDebugMacro(<<" Distribution invisible, because numColors == this->myNbValues.size()"); - if ( distrVisibility && GetDistributionVisibility() ) { + if ( distrVisibility && GetDistributionVisibility() ) + { for ( i = 0 ; i < (int)myNbValues.size(); i++ ) { if ( myNbValues[i] ) { numPositiveVal++; @@ -349,16 +352,21 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) this->myDistribution->SetPolys(distrPolys); distrPts->Delete(); distrPolys->Delete(); - if ( myDistributionColoringType == SMESH_MULTICOLOR_TYPE ) { + if ( myDistributionColoringType == SMESH_MULTICOLOR_TYPE ) + { distColors = vtkUnsignedCharArray::New(); distColors->SetNumberOfComponents(3); distColors->SetNumberOfTuples(numPositiveVal); this->myDistribution->GetCellData()->SetScalars(distColors); distColors->Delete(); - } else if( myDistributionColoringType == SMESH_MONOCOLOR_TYPE ){ + } + else if( myDistributionColoringType == SMESH_MONOCOLOR_TYPE ) + { this->myDistribution->GetCellData()->SetScalars(NULL); } - } else { + } + else + { myDistribution->Reset(); } // rnv end @@ -373,22 +381,22 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) // get the viewport size in display coordinates int *barOrigin, barWidth, barHeight, distrHeight; barOrigin = this->PositionCoordinate->GetComputedViewportValue(viewport); - size[0] = + size[0] = this->Position2Coordinate->GetComputedViewportValue(viewport)[0] - barOrigin[0]; - size[1] = + size[1] = this->Position2Coordinate->GetComputedViewportValue(viewport)[1] - barOrigin[1]; this->LastOrigin[0] = barOrigin[0]; this->LastOrigin[1] = barOrigin[1]; this->LastSize[0] = size[0]; this->LastSize[1] = size[1]; - + // Update all the composing objects this->TitleActor->SetProperty(this->GetProperty()); this->TitleMapper->SetInput(this->Title); if (this->TitleTextProperty->GetMTime() > this->BuildTime) - { + { // Shallow copy here so that the size of the title prop is not affected // by the automatic adjustment of its text mapper's size (i.e. its // mapper's text property is identical except for the font size @@ -397,17 +405,17 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) // the title and label text prop to be the same. this->TitleMapper->GetTextProperty()->ShallowCopy(this->TitleTextProperty); this->TitleMapper->GetTextProperty()->SetJustificationToCentered(); - } - + } + // find the best size for the title font int titleSize[2]; this->SizeTitle(titleSize, size, viewport); - + // find the best size for the ticks int labelSize[2]; this->AllocateAndSizeLabels(labelSize, size, viewport,range); this->NumberOfLabelsBuilt = this->NumberOfLabels; - + // generate points double x[3]; x[2] = 0.0; double delta, itemH, shrink; @@ -426,7 +434,7 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) barHeight = (int)(0.86*size[1]); delta=(double)barHeight/numColors; - + for ( i=0; iSetPoint(2*i+1,x); } - if(GetDistributionVisibility() && distrVisibility) { - // Distribution points + if ( GetDistributionVisibility() && distrVisibility ) { + // Distribution points shrink = delta*SHRINK_COEF; vtkIdType distPtsId=0; vtkIdType distPtsIds[4]; - for(i=0; iSetPoint(distPtsId++,x); @@ -462,7 +470,7 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) x[1] = i*delta+delta-shrink; // third point of polygon (quadrangle) - x[0] = 0; + x[0] = 0; distPtsIds[3] = distPtsId; distrPts->SetPoint(distPtsId++,x); @@ -475,19 +483,19 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) distrPolys->InsertNextCell(4,distPtsIds); } } - } + } } // rnv end else { barWidth = size[0]; - + // rnv begin // Customization of the vtkScalarBarActor to show distribution histogram. double coef1, delimeter=0.0; - if(GetDistributionVisibility() && distrVisibility) { + if ( GetDistributionVisibility() && distrVisibility ) { coef1=0.62; distrHeight = (int)((coef1/2)*size[1]); - //delimeter between distribution diagram and scalar bar + //delimeter between distribution diagram and scalar bar delimeter=0.02*size[1]; } else { @@ -495,64 +503,64 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) barHeight = (int)(coef1*size[1]); distrHeight = 0; } - + barHeight = (int)(coef1*size[1]); - + delta=(double)barWidth/numColors; - for (i=0; iSetPoint(2*i,x); + pts->SetPoint(2*i,x); x[1] = distrHeight + delimeter; pts->SetPoint(2*i+1,x); } - - if(GetDistributionVisibility() && distrVisibility) { - // Distribution points + + if ( GetDistributionVisibility() && distrVisibility ) { + // Distribution points shrink = delta*SHRINK_COEF; vtkIdType distPtsId=0; vtkIdType distPtsIds[4]; - for(i=0; iSetPoint(distPtsId++,x); - + // second point of polygon (quadrangle) - x[0] = i*delta+shrink; + x[0] = i*delta+shrink; x[1] = itemH; distPtsIds[3] = distPtsId; distrPts->SetPoint(distPtsId++,x); - + // third point of polygon (quadrangle) - x[0] = i*delta+delta-shrink; + x[0] = i*delta+delta-shrink; x[1] = 0; distPtsIds[1] = distPtsId; distrPts->SetPoint(distPtsId++,x); - + // fourth point of polygon (quadrangle) - x[0] = i*delta+delta-shrink; + x[0] = i*delta+delta-shrink; x[1] = itemH; distPtsIds[2] = distPtsId; distrPts->SetPoint(distPtsId++,x); - + // Add polygon into poly data distrPolys->InsertNextCell(4,distPtsIds); } - } + } } // rnv end } - + //polygons & cell colors unsigned char *rgba, *rgb; vtkIdType ptIds[4], dcCount=0; - for (i=0; iInsertNextCell(4,ptIds); if ( isLogTable ) - { - double rgbval = log10(range[0]) + + { + double rgbval = log10(range[0]) + i*(log10(range[1])-log10(range[0]))/(numColors -1); rgba = lut->MapValue(pow(10.0,rgbval)); - } + } else - { + { rgba = lut->MapValue(range[0] + (range[1] - range[0])* ((double)i /(numColors-1.0))); - } + } rgb = colors->GetPointer(3*i); //write into array directly rgb[0] = rgba[0]; rgb[1] = rgba[1]; rgb[2] = rgba[2]; - + // rnv begin // Customization of the vtkScalarBarActor to show distribution histogram. if ( myDistributionColoringType == SMESH_MULTICOLOR_TYPE && @@ -594,62 +602,62 @@ int SMESH_ScalarBarActor::RenderOpaqueGeometry(vtkViewport *viewport) // Now position everything properly // double val; - if (this->Orientation == VTK_ORIENT_VERTICAL) - { + if ( this->Orientation == VTK_ORIENT_VERTICAL ) + { int sizeTextData[2]; - + // center the title this->TitleActor->SetPosition(size[0]/2, 0.9*size[1]); - - for (i=0; i < this->NumberOfLabels; i++) + + for ( i = 0; i < this->NumberOfLabels; i++ ) + { + if ( this->NumberOfLabels > 1 ) { - if (this->NumberOfLabels > 1) - { val = (double)i/(this->NumberOfLabels-1) *barHeight; - } - else - { + } + else + { val = 0.5*barHeight; - } + } this->TextMappers[i]->GetSize(viewport,sizeTextData); this->TextMappers[i]->GetTextProperty()->SetJustificationToLeft(); this->TextActors[i]->SetPosition(barWidth+3, val - sizeTextData[1]/2); - } } + } else - { - this->TitleActor->SetPosition(size[0]/2, + { + this->TitleActor->SetPosition(size[0]/2, barHeight + labelSize[1] + 0.1*size[1]); - for (i=0; i < this->NumberOfLabels; i++) - { + for ( i = 0; i < this->NumberOfLabels; i++ ) + { this->TextMappers[i]->GetTextProperty()->SetJustificationToCentered(); if (this->NumberOfLabels > 1) - { + { val = (double)i/(this->NumberOfLabels-1) * barWidth; - } + } else - { + { val = 0.5*barWidth; - } - this->TextActors[i]->SetPosition(val, barHeight + 0.05*size[1]); } + this->TextActors[i]->SetPosition(val, barHeight + 0.05*size[1]); } + } this->BuildTime.Modified(); - } + } // Everything is built, just have to render - if (this->Title != NULL) - { + if ( this->Title != NULL ) + { renderedSomething += this->TitleActor->RenderOpaqueGeometry(viewport); - } + } this->ScalarBarActor->RenderOpaqueGeometry(viewport); this->myDistributionActor->RenderOpaqueGeometry(viewport); - for (i=0; iNumberOfLabels; i++) - { + for ( i = 0; i < this->NumberOfLabels; i++ ) + { renderedSomething += this->TextActors[i]->RenderOpaqueGeometry(viewport); - } + } renderedSomething = (renderedSomething > 0)?(1):(0); @@ -662,50 +670,50 @@ void SMESH_ScalarBarActor::PrintSelf(ostream& os, vtkIndent indent) this->Superclass::PrintSelf(os,indent); if ( this->LookupTable ) - { + { os << indent << "Lookup Table:\n"; this->LookupTable->PrintSelf(os,indent.GetNextIndent()); - } + } else - { + { os << indent << "Lookup Table: (none)\n"; - } + } if (this->TitleTextProperty) - { + { os << indent << "Title Text Property:\n"; this->TitleTextProperty->PrintSelf(os,indent.GetNextIndent()); - } + } else - { + { os << indent << "Title Text Property: (none)\n"; - } + } if (this->LabelTextProperty) - { + { os << indent << "Label Text Property:\n"; this->LabelTextProperty->PrintSelf(os,indent.GetNextIndent()); - } + } else - { + { os << indent << "Label Text Property: (none)\n"; - } + } os << indent << "Title: " << (this->Title ? this->Title : "(none)") << "\n"; - os << indent << "Maximum Number Of Colors: " + os << indent << "Maximum Number Of Colors: " << this->MaximumNumberOfColors << "\n"; os << indent << "Number Of Labels: " << this->NumberOfLabels << "\n"; os << indent << "Number Of Labels Built: " << this->NumberOfLabelsBuilt << "\n"; os << indent << "Orientation: "; if ( this->Orientation == VTK_ORIENT_HORIZONTAL ) - { + { os << "Horizontal\n"; - } + } else - { + { os << "Vertical\n"; - } + } os << indent << "Label Format: " << this->LabelFormat << "\n"; } @@ -715,7 +723,7 @@ void SMESH_ScalarBarActor::ShallowCopy(vtkProp *prop) { SMESH_ScalarBarActor *a = SMESH_ScalarBarActor::SafeDownCast(prop); if ( a != NULL ) - { + { this->SetPosition2(a->GetPosition2()); this->SetLookupTable(a->GetLookupTable()); this->SetMaximumNumberOfColors(a->GetMaximumNumberOfColors()); @@ -724,25 +732,25 @@ void SMESH_ScalarBarActor::ShallowCopy(vtkProp *prop) this->SetTitleTextProperty(a->GetTitleTextProperty()); this->SetLabelFormat(a->GetLabelFormat()); this->SetTitle(a->GetTitle()); - this->GetPositionCoordinate()->SetCoordinateSystem( - a->GetPositionCoordinate()->GetCoordinateSystem()); - this->GetPositionCoordinate()->SetValue( - a->GetPositionCoordinate()->GetValue()); - this->GetPosition2Coordinate()->SetCoordinateSystem( - a->GetPosition2Coordinate()->GetCoordinateSystem()); - this->GetPosition2Coordinate()->SetValue( - a->GetPosition2Coordinate()->GetValue()); - } + this->GetPositionCoordinate()->SetCoordinateSystem + (a->GetPositionCoordinate()->GetCoordinateSystem()); + this->GetPositionCoordinate()->SetValue + (a->GetPositionCoordinate()->GetValue()); + this->GetPosition2Coordinate()->SetCoordinateSystem + (a->GetPosition2Coordinate()->GetCoordinateSystem()); + this->GetPosition2Coordinate()->SetValue + (a->GetPosition2Coordinate()->GetValue()); + } // Now do superclass this->vtkActor2D::ShallowCopy(prop); } //---------------------------------------------------------------------------- -void SMESH_ScalarBarActor::AllocateAndSizeLabels(int *labelSize, - int *size, - vtkViewport *viewport, - double *range) +void SMESH_ScalarBarActor::AllocateAndSizeLabels(int *labelSize, + int *size, + vtkViewport *viewport, + double *range) { labelSize[0] = labelSize[1] = 0; @@ -753,55 +761,55 @@ void SMESH_ScalarBarActor::AllocateAndSizeLabels(int *labelSize, double val; int i; - + // TODO: this should be optimized, maybe by keeping a list of // allocated mappers, in order to avoid creation/destruction of // their underlying text properties (i.e. each time a mapper is // created, text properties are created and shallow-assigned a font size // which value might be "far" from the target font size). - // is this a vtkLookupTable or a subclass of vtkLookupTable + // is this a vtkLookupTable or a subclass of vtkLookupTable // with its scale set to log vtkLookupTable *LUT = vtkLookupTable::SafeDownCast( this->LookupTable ); int isLogTable = 0; if ( LUT ) - { + { if ( LUT->GetScale() == VTK_SCALE_LOG10 ) - { - isLogTable = 1; - } + { + isLogTable = 1; } + } - for (i=0; i < this->NumberOfLabels; i++) - { + for ( i = 0; i < this->NumberOfLabels; i++ ) + { this->TextMappers[i] = vtkTextMapper::New(); if ( isLogTable ) - { + { double lval; - if (this->NumberOfLabels > 1) - { + if ( this->NumberOfLabels > 1 ) + { lval = log10(range[0]) + (double)i/(this->NumberOfLabels-1) * (log10(range[1])-log10(range[0])); - } + } else - { + { lval = log10(range[0]) + 0.5*(log10(range[1])-log10(range[0])); - } - val = pow(10.0,lval); } + val = pow(10.0,lval); + } else + { + if ( this->NumberOfLabels > 1 ) { - if (this->NumberOfLabels > 1) - { - val = range[0] + + val = range[0] + (double)i/(this->NumberOfLabels-1) * (range[1]-range[0]); - } + } else - { + { val = range[0] + 0.5*(range[1]-range[0]); - } } + } sprintf(string, this->LabelFormat, val); this->TextMappers[i]->SetInput(string); @@ -812,64 +820,63 @@ void SMESH_ScalarBarActor::AllocateAndSizeLabels(int *labelSize, // which will be modified later). This allows text actors to // share the same text property, and in that case specifically allows // the title and label text prop to be the same. - this->TextMappers[i]->GetTextProperty()->ShallowCopy( - this->LabelTextProperty); + this->TextMappers[i]->GetTextProperty()->ShallowCopy(this->LabelTextProperty); this->TextActors[i] = vtkActor2D::New(); this->TextActors[i]->SetMapper(this->TextMappers[i]); this->TextActors[i]->SetProperty(this->GetProperty()); this->TextActors[i]->GetPositionCoordinate()-> SetReferenceCoordinate(this->PositionCoordinate); - } + } - if (this->NumberOfLabels) - { + if ( this->NumberOfLabels ) + { int targetWidth, targetHeight; // rnv begin // Customization of the vtkScalarBarActor to show distribution histogram. bool distrVisibility = ( this->MaximumNumberOfColors == (int) this->myNbValues.size() ); double coef; - if( GetDistributionVisibility() && distrVisibility ) - if(this->Orientation == VTK_ORIENT_VERTICAL) + if ( GetDistributionVisibility() && distrVisibility ) + if ( this->Orientation == VTK_ORIENT_VERTICAL ) coef = 0.4; - else + else coef = 0.18; - else - if(this->Orientation == VTK_ORIENT_VERTICAL) + else + if (this->Orientation == VTK_ORIENT_VERTICAL ) coef = 0.6; - else + else coef=0.25; if ( this->Orientation == VTK_ORIENT_VERTICAL ) - { + { targetWidth = (int)(coef*size[0]); targetHeight = (int)(0.86*size[1]/this->NumberOfLabels); - } + } else - { + { targetWidth = (int)(size[0]*0.8/this->NumberOfLabels); targetHeight = (int)(coef*size[1]); - } - // rnv end - - vtkTextMapper::SetMultipleConstrainedFontSize(viewport, - targetWidth, - targetHeight, - this->TextMappers, - this->NumberOfLabels, - labelSize); } + // rnv end + + vtkTextMapper::SetMultipleConstrainedFontSize( viewport, + targetWidth, + targetHeight, + this->TextMappers, + this->NumberOfLabels, + labelSize ); + } } //---------------------------------------------------------------------------- -void SMESH_ScalarBarActor::SizeTitle(int *titleSize, - int *size, +void SMESH_ScalarBarActor::SizeTitle(int *titleSize, + int *size, vtkViewport *viewport) { titleSize[0] = titleSize[1] = 0; - if (this->Title == NULL || !strlen(this->Title)) + if ( this->Title == NULL || !strlen(this->Title) ) { return; } @@ -902,36 +909,43 @@ void SMESH_ScalarBarActor::SizeTitle(int *titleSize, /*--------------------------------------------------------------------------*/ -void SMESH_ScalarBarActor::SetDistributionVisibility(int flag) { - myDistributionActor->SetVisibility(flag); +void SMESH_ScalarBarActor::SetDistributionVisibility( int flag ) +{ + myDistributionActor->SetVisibility( flag ); Modified(); } /*--------------------------------------------------------------------------*/ -int SMESH_ScalarBarActor::GetDistributionVisibility() { +int SMESH_ScalarBarActor::GetDistributionVisibility() +{ return myDistributionActor->GetVisibility(); } -void SMESH_ScalarBarActor::SetDistribution(std::vector theNbValues) { +void SMESH_ScalarBarActor::SetDistribution( const std::vector& theNbValues ) +{ myNbValues = theNbValues; -} +} -void SMESH_ScalarBarActor::SetDistributionColor (double rgb[3]) { +void SMESH_ScalarBarActor::SetDistributionColor( double rgb[3] ) +{ myDistributionActor->GetProperty()->SetColor(rgb); Modified(); } -void SMESH_ScalarBarActor::GetDistributionColor (double rgb[3]) { +void SMESH_ScalarBarActor::GetDistributionColor( double rgb[3] ) +{ myDistributionActor->GetProperty()->GetColor(rgb); } -void SMESH_ScalarBarActor::SetTitleOnlyVisibility( bool theTitleOnlyVisibility) { +void SMESH_ScalarBarActor::SetTitleOnlyVisibility( bool theTitleOnlyVisibility ) +{ myTitleOnlyVisibility = theTitleOnlyVisibility; } -bool SMESH_ScalarBarActor::GetTitleOnlyVisibility() { +bool SMESH_ScalarBarActor::GetTitleOnlyVisibility() +{ return myTitleOnlyVisibility; } diff --git a/src/OBJECT/SMESH_ScalarBarActor.h b/src/OBJECT/SMESH_ScalarBarActor.h index ffc91bb11..c5e721b23 100644 --- a/src/OBJECT/SMESH_ScalarBarActor.h +++ b/src/OBJECT/SMESH_ScalarBarActor.h @@ -174,7 +174,7 @@ class SMESHOBJECT_EXPORT SMESH_ScalarBarActor: public vtkActor2D { virtual int GetDistributionVisibility(); // Description: // Set distribution - virtual void SetDistribution(std::vector theNbValues); + virtual void SetDistribution(const std::vector& theNbValues); // Description: // Set distribution coloring type (SMESH_MONOCOLOR_TYPE or SMESH_MULTICOLOR_TYPE) diff --git a/src/SMDS/SMDS_Mesh.cxx b/src/SMDS/SMDS_Mesh.cxx index 5d8294700..7c8e66327 100644 --- a/src/SMDS/SMDS_Mesh.cxx +++ b/src/SMDS/SMDS_Mesh.cxx @@ -2458,6 +2458,49 @@ const SMDS_MeshElement* SMDS_Mesh::FindElement (const vector& nodes, + std::vector& foundElems, + const SMDSAbs_ElementType type) +{ + // chose a node with minimal number of inverse elements + const SMDS_MeshNode* n0 = nodes[0]; + int minNbInverse = n0 ? n0->NbInverseElements( type ) : 1000; + for ( size_t i = 1; i < nodes.size(); ++i ) + if ( nodes[i] && nodes[i]->NbInverseElements( type ) < minNbInverse ) + { + n0 = nodes[i]; + minNbInverse = n0->NbInverseElements( type ); + } + + foundElems.clear(); + if ( n0 ) + { + foundElems.reserve( minNbInverse ); + SMDS_ElemIteratorPtr eIt = n0->GetInverseElementIterator( type ); + while ( eIt->more() ) + { + const SMDS_MeshElement* e = eIt->next(); + bool includeAll = true; + for ( size_t i = 0; i < nodes.size() && includeAll; ++i ) + if ( nodes[i] != n0 && e->GetNodeIndex( nodes[i] ) < 0 ) + includeAll = false; + if ( includeAll ) + foundElems.push_back( e ); + } + } + return foundElems.size(); +} + //======================================================================= //function : DumpNodes //purpose : diff --git a/src/SMDS/SMDS_Mesh.hxx b/src/SMDS/SMDS_Mesh.hxx index b1af99a1d..d1d0b1120 100644 --- a/src/SMDS/SMDS_Mesh.hxx +++ b/src/SMDS/SMDS_Mesh.hxx @@ -689,6 +689,9 @@ public: static const SMDS_MeshElement* FindElement(const std::vector& nodes, const SMDSAbs_ElementType type=SMDSAbs_All, const bool noMedium=true); + static int GetElementsByNodes(const std::vector& nodes, + std::vector& foundElems, + const SMDSAbs_ElementType type=SMDSAbs_All); /*! * \brief Raise an exception if free memory (ram+swap) too low diff --git a/src/SMESH/SMESH_Algo.cxx b/src/SMESH/SMESH_Algo.cxx index 035c23c59..4e4002ee3 100644 --- a/src/SMESH/SMESH_Algo.cxx +++ b/src/SMESH/SMESH_Algo.cxx @@ -553,28 +553,35 @@ bool SMESH_Algo::IsStraight( const TopoDS_Edge & E, case GeomAbs_Hyperbola: case GeomAbs_Parabola: return false; - // case GeomAbs_BezierCurve: - // case GeomAbs_BSplineCurve: - // case GeomAbs_OtherCurve: + // case GeomAbs_BezierCurve: + // case GeomAbs_BSplineCurve: + // case GeomAbs_OtherCurve: default:; } - const double f = curve.FirstParameter(); - const double l = curve.LastParameter(); - const gp_Pnt pf = curve.Value( f ); - const gp_Pnt pl = curve.Value( l ); - const gp_Vec v1( pf, pl ); - const double v1Len = v1.Magnitude(); - if ( v1Len < std::numeric_limits< double >::min() ) + + // evaluate how far from a straight line connecting the curve ends + // stand internal points of the curve + double f = curve.FirstParameter(); + double l = curve.LastParameter(); + gp_Pnt pf = curve.Value( f ); + gp_Pnt pl = curve.Value( l ); + gp_Vec lineVec( pf, pl ); + double lineLen2 = lineVec.SquareMagnitude(); + if ( lineLen2 < std::numeric_limits< double >::min() ) return false; // E seems closed - const double tol = Min( 10 * curve.Tolerance(), v1Len * 1e-2 ); + + double edgeTol = 10 * curve.Tolerance(); + double lenTol2 = lineLen2 * 1e-4; + double tol2 = Min( edgeTol * edgeTol, lenTol2 ); + const double nbSamples = 7; for ( int i = 0; i < nbSamples; ++i ) { - const double r = ( i + 1 ) / nbSamples; - const gp_Pnt pi = curve.Value( f * r + l * ( 1 - r )); - const gp_Vec vi( pf, pi ); - const double h = 0.5 * v1.Crossed( vi ).Magnitude() / v1Len; - if ( h > tol ) + double r = ( i + 1 ) / nbSamples; + gp_Pnt pi = curve.Value( f * r + l * ( 1 - r )); + gp_Vec vi( pf, pi ); + double h2 = lineVec.Crossed( vi ).SquareMagnitude() / lineLen2; + if ( h2 > tol2 ) return false; } return true; @@ -852,6 +859,16 @@ bool SMESH_Algo::Compute(SMESH_Mesh & /*aMesh*/, SMESH_MesherHelper* /*aHelper*/ return error( COMPERR_BAD_INPUT_MESH, "Mesh built on shape expected"); } +//======================================================================= +//function : IsApplicableToShape +//purpose : Return true if the algorithm can mesh a given shape +//======================================================================= + +bool SMESH_Algo::IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const +{ + return true; +} + //======================================================================= //function : CancelCompute //purpose : Sets _computeCanceled to true. It's usage depends on diff --git a/src/SMESH/SMESH_Algo.hxx b/src/SMESH/SMESH_Algo.hxx index d4c8f3e87..ba8c2d356 100644 --- a/src/SMESH/SMESH_Algo.hxx +++ b/src/SMESH/SMESH_Algo.hxx @@ -165,6 +165,15 @@ class SMESH_EXPORT SMESH_Algo : public SMESH_Hypothesis */ virtual bool Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelper); + /*! + * \brief Return true if the algorithm can mesh a given shape + * \param [in] aShape - shape to check + * \param [in] toCheckAll - if true, this check returns OK if all shapes are OK, + * else, returns OK if at least one shape is OK + * \retval bool - \c true by default + */ + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const; + /*! * \brief Sets _computeCanceled to true. It's usage depends on * implementation of a particular mesher. diff --git a/src/SMESH/SMESH_MeshEditor.cxx b/src/SMESH/SMESH_MeshEditor.cxx index 1f343b1c4..c477e729e 100644 --- a/src/SMESH/SMESH_MeshEditor.cxx +++ b/src/SMESH/SMESH_MeshEditor.cxx @@ -7482,8 +7482,8 @@ void SMESH_MeshEditor::MergeNodes (TListOfListOfNodes & theGroupsOfNodes, } // Remove bad elements, then equal nodes (order important) - Remove( rmElemIds, false ); - Remove( rmNodeIds, true ); + Remove( rmElemIds, /*isNodes=*/false ); + Remove( rmNodeIds, /*isNodes=*/true ); return; } @@ -7553,14 +7553,14 @@ bool SMESH_MeshEditor::applyMerge( const SMDS_MeshElement* elem, toRemove = true; nbResElems = 0; - if ( elem->IsQuadratic() && newElemDefs[0].myType == SMDSAbs_Face && nbNodes > 6 ) + if ( newElemDefs[0].myIsQuad && newElemDefs[0].myType == SMDSAbs_Face && nbNodes > 6 ) { // if corner nodes stick, remove medium nodes between them from uniqueNodes int nbCorners = nbNodes / 2; for ( int iCur = 0; iCur < nbCorners; ++iCur ) { - int iPrev = ( iCur + 1 ) % nbCorners; - if ( curNodes[ iCur ] == curNodes[ iPrev ] ) // corners stick + int iNext = ( iCur + 1 ) % nbCorners; + if ( curNodes[ iCur ] == curNodes[ iNext ] ) // corners stick { int iMedium = iCur + nbCorners; vector< const SMDS_MeshNode* >::iterator i = @@ -7711,11 +7711,9 @@ bool SMESH_MeshEditor::applyMerge( const SMDS_MeshElement* elem, // | | // +---+---+ // 0 7 3 - if (( nbUniqueNodes == 7 && nbRepl == 2 && iRepl[1] != 8 ) && - (( iRepl[0] == 1 && iRepl[1] == 4 && curNodes[1] == curNodes[0] ) || - ( iRepl[0] == 2 && iRepl[1] == 5 && curNodes[2] == curNodes[1] ) || - ( iRepl[0] == 3 && iRepl[1] == 6 && curNodes[3] == curNodes[2] ) || - ( iRepl[0] == 3 && iRepl[1] == 7 && curNodes[3] == curNodes[0] ))) + if ( nbUniqueNodes == 7 && + iRepl[0] < 4 && + ( nbRepl == 1 || iRepl[1] != 8 )) { toRemove = false; } diff --git a/src/SMESH/SMESH_MesherHelper.cxx b/src/SMESH/SMESH_MesherHelper.cxx index e03128dc8..004c856f1 100644 --- a/src/SMESH/SMESH_MesherHelper.cxx +++ b/src/SMESH/SMESH_MesherHelper.cxx @@ -2669,7 +2669,7 @@ bool SMESH_MesherHelper::LoadNodeColumns(TParam2ColumnMap & theParam2 } // nb rows of nodes - size_t prevNbRows = theParam2ColumnMap.begin()->second.size(); // current, at least 1 here + size_t prevNbRows = theParam2ColumnMap.begin()->second.size(); // current, at least 1 here size_t expectNbRows = faceSubMesh->NbElements() / ( theParam2ColumnMap.size()-1 ); // to be added // fill theParam2ColumnMap column by column by passing from nodes on diff --git a/src/SMESH/SMESH_Pattern.cxx b/src/SMESH/SMESH_Pattern.cxx index 0f14f5ef6..db6e60e8c 100644 --- a/src/SMESH/SMESH_Pattern.cxx +++ b/src/SMESH/SMESH_Pattern.cxx @@ -602,8 +602,7 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, Extrema_GenExtPS projector; GeomAdaptor_Surface aSurface( BRep_Tool::Surface( face )); - if ( theProject || needProject ) - projector.Initialize( aSurface, 20,20, 1e-5,1e-5 ); + projector.Initialize( aSurface, 20,20, 1e-5,1e-5 ); int iPoint = 0; TNodePointIDMap nodePointIDMap; @@ -690,8 +689,28 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, myPoints.resize( nbNodes ); + // care of INTERNAL VERTEXes + TopExp_Explorer vExp( face, TopAbs_VERTEX, TopAbs_EDGE ); + for ( ; vExp.More(); vExp.Next() ) + { + const SMDS_MeshNode* node = + SMESH_Algo::VertexNode( TopoDS::Vertex( vExp.Current()), aMeshDS ); + if ( !node || node->NbInverseElements( SMDSAbs_Face ) == 0 ) + continue; + myPoints.resize( ++nbNodes ); + list< TPoint* > & fPoints = getShapePoints( face ); + nodePointIDMap.insert( make_pair( node, iPoint )); + TPoint* p = &myPoints[ iPoint++ ]; + fPoints.push_back( p ); + gp_XY uv = helper.GetNodeUV( face, node ); + p->myInitUV.SetCoord( uv.X(), uv.Y() ); + p->myInitXYZ.SetCoord( p->myInitUV.X(), p->myInitUV.Y(), 0 ); + } + // Load U of points on edges + Bnd_Box2d edgesUVBox; + list::iterator nbEinW = myNbKeyPntInBoundary.begin(); int iE = 0; vector< TopoDS_Edge > eVec; @@ -762,6 +781,7 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, else keyPoint->myInitUV = C2d->Value( isForward ? f : l ).XY(); keyPoint->myInitXYZ.SetCoord (keyPoint->myInitUV.X(), keyPoint->myInitUV.Y(), 0); + edgesUVBox.Add( gp_Pnt2d( keyPoint->myInitUV )); } } if ( !vPoint->empty() ) @@ -841,6 +861,7 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, p->myInitUV = C2d->Value( u ).XY(); } p->myInitXYZ.SetCoord( p->myInitUV.X(), p->myInitUV.Y(), 0 ); + edgesUVBox.Add( gp_Pnt2d( p->myInitUV )); unIt++; unRIt++; iPoint++; } @@ -866,6 +887,7 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, else keyPoint->myInitUV = C2d->Value( isForward ? l : f ).XY(); keyPoint->myInitXYZ.SetCoord( keyPoint->myInitUV.X(), keyPoint->myInitUV.Y(), 0 ); + edgesUVBox.Add( gp_Pnt2d( keyPoint->myInitUV )); } } if ( !vPoint->empty() ) @@ -906,7 +928,7 @@ bool SMESH_Pattern::Load (SMESH_Mesh* theMesh, nodePointIDMap.insert( make_pair( node, iPoint )); TPoint* p = &myPoints[ iPoint++ ]; fPoints.push_back( p ); - if ( theProject ) + if ( theProject || edgesUVBox.IsOut( p->myInitUV ) ) p->myInitUV = project( node, projector ); else { const SMDS_FacePosition* pos = diff --git a/src/SMESHGUI/SMESHGUI.cxx b/src/SMESHGUI/SMESHGUI.cxx index 526b6c650..1145efb17 100644 --- a/src/SMESHGUI/SMESHGUI.cxx +++ b/src/SMESHGUI/SMESHGUI.cxx @@ -630,7 +630,7 @@ namespace }; // is typeMsg complete? (compilation failure mains that enum SMDSAbs_EntityType changed) const int nbTypes = sizeof( typeMsg ) / sizeof( const char* ); - int _assert[( nbTypes == SMESH::Entity_Last ) ? 2 : -1 ]; _assert[0]=_assert[1]; + int _assert[( nbTypes == SMESH::Entity_Last ) ? 2 : -1 ]; _assert[0]=_assert[1]=0; QString andStr = " " + QObject::tr("SMESH_AND") + " ", comma(", "); for ( size_t iType = 0; iType < presentNotSupported.size(); ++iType ) { diff --git a/src/SMESHGUI/SMESHGUI_DisplayEntitiesDlg.cxx b/src/SMESHGUI/SMESHGUI_DisplayEntitiesDlg.cxx index 5af7f7815..14088fe1a 100644 --- a/src/SMESHGUI/SMESHGUI_DisplayEntitiesDlg.cxx +++ b/src/SMESHGUI/SMESHGUI_DisplayEntitiesDlg.cxx @@ -97,7 +97,7 @@ SMESHGUI_DisplayEntitiesDlg::SMESHGUI_DisplayEntitiesDlg( QWidget* parent ) hl->addWidget( nb0DElemsLab, 0, 1 ); my0DElemsTB->setEnabled( nbElements ); nb0DElemsLab->setEnabled( nbElements ); - myNbTypes += ( nbElements > 0 ); + myNbTypes = ( nbElements > 0 ); // Edges nbElements = myActor ? myActor->GetObject()->GetNbEntities( SMDSAbs_Edge ) : aMesh->NbEdges(); diff --git a/src/SMESHGUI/SMESHGUI_MeshInfo.cxx b/src/SMESHGUI/SMESHGUI_MeshInfo.cxx index c2932c203..299ba97c0 100644 --- a/src/SMESHGUI/SMESHGUI_MeshInfo.cxx +++ b/src/SMESHGUI/SMESHGUI_MeshInfo.cxx @@ -77,22 +77,24 @@ namespace { -const int SPACING = 6; -const int MARGIN = 9; -const int MAXITEMS = 10; -const int GROUPS_ID = 100; -const int SUBMESHES_ID = 200; -const int SPACING_INFO = 2; - -enum InfoRole { - TypeRole = Qt::UserRole + 10, - IdRole, -}; - -enum InfoType { - NodeConnectivity = 100, - ElemConnectivity, -}; + const int SPACING = 6; + const int MARGIN = 9; + const int MAXITEMS = 10; + const int GROUPS_ID = 100; + const int SUBMESHES_ID = 200; + const int SPACING_INFO = 2; + + const char* id_preview_resource = "id_preview_resource"; + + enum InfoRole { + TypeRole = Qt::UserRole + 10, + IdRole, + }; + + enum InfoType { + NodeConnectivity = 100, + ElemConnectivity, + }; } // namesapce /*! @@ -2942,6 +2944,8 @@ SMESHGUI_MeshInfoDlg::SMESHGUI_MeshInfoDlg( QWidget* parent, int page ) connect( myElemInfo, SIGNAL( itemInfo( int )), this, SLOT( showItemInfo( int ))); connect( myElemInfo, SIGNAL( itemInfo( QString )), this, SLOT( showItemInfo( QString ))); + myIDPreviewCheck->setChecked( SMESHGUI::resourceMgr()->booleanValue( "SMESH", id_preview_resource, false )); + updateSelection(); } @@ -3189,6 +3193,7 @@ void SMESHGUI_MeshInfoDlg::idChanged() void SMESHGUI_MeshInfoDlg::idPreviewChange( bool isOn ) { myIDPreview->SetPointsLabeled( isOn && !myID->text().simplified().isEmpty() ); + SMESHGUI::resourceMgr()->setValue("SMESH", id_preview_resource, isOn ); if ( SVTK_ViewWindow* aViewWindow = SMESH::GetViewWindow() ) aViewWindow->Repaint(); } diff --git a/src/SMESHUtils/CMakeLists.txt b/src/SMESHUtils/CMakeLists.txt index 3641d913c..c84715eed 100644 --- a/src/SMESHUtils/CMakeLists.txt +++ b/src/SMESHUtils/CMakeLists.txt @@ -67,6 +67,7 @@ SET(SMESHUtils_HEADERS SMESH_MeshAlgos.hxx SMESH_MAT2d.hxx SMESH_ControlPnt.hxx + SMESH_Delaunay.hxx ) # --- sources --- @@ -84,6 +85,7 @@ SET(SMESHUtils_SOURCES SMESH_FreeBorders.cxx SMESH_ControlPnt.cxx SMESH_DeMerge.cxx + SMESH_Delaunay.cxx ) # --- rules --- diff --git a/src/SMESHUtils/SMESH_Delaunay.cxx b/src/SMESHUtils/SMESH_Delaunay.cxx new file mode 100644 index 000000000..e53bfd8fc --- /dev/null +++ b/src/SMESHUtils/SMESH_Delaunay.cxx @@ -0,0 +1,373 @@ +// Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// File : SMESH_Delaunay.cxx +// Created : Wed Apr 19 15:41:15 2017 +// Author : Edward AGAPOV (eap) + +#include "SMESH_Delaunay.hxx" + +#include "SMESH_Comment.hxx" +#include "SMESH_File.hxx" +#include "SMESH_MeshAlgos.hxx" + +#include +#include + +//================================================================================ +/*! + * \brief Construct a Delaunay triangulation of given boundary nodes + * \param [in] boundaryNodes - vector of nodes of a wire + * \param [in] face - the face + * \param [in] faceID - the face ID + * \param [in] nbNodesToVisit - nb of non-marked nodes on the face + */ +//================================================================================ + +SMESH_Delaunay::SMESH_Delaunay(const std::vector< const UVPtStructVec* > & boundaryNodes, + const TopoDS_Face& face, + const int faceID) + : _face( face ), _faceID( faceID ), _scale( 1., 1. ) +{ + // compute _scale + BRepAdaptor_Surface surf( face ); + if ( surf.GetType() != GeomAbs_Plane ) + { + const int nbDiv = 100; + const double uRange = surf.LastUParameter() - surf.FirstUParameter(); + const double vRange = surf.LastVParameter() - surf.FirstVParameter(); + const double dU = uRange / nbDiv; + const double dV = vRange / nbDiv; + double u = surf.FirstUParameter(), v = surf.FirstVParameter(); + gp_Pnt p0U = surf.Value( u, v ), p0V = p0U; + double lenU = 0, lenV = 0; + for ( ; u < surf.LastUParameter(); u += dU, v += dV ) + { + gp_Pnt p1U = surf.Value( u, surf.FirstVParameter() ); + lenU += p1U.Distance( p0U ); + p0U = p1U; + gp_Pnt p1V = surf.Value( surf.FirstUParameter(), v ); + lenV += p1V.Distance( p0V ); + p0V = p1V; + } + _scale.SetCoord( lenU / uRange, lenV / vRange ); + } + + // count boundary points + int iP = 1, nbP = 0; + for ( size_t iW = 0; iW < boundaryNodes.size(); ++iW ) // loop on wires + { + nbP += boundaryNodes[iW]->size(); + if ( boundaryNodes[iW]->front().node == boundaryNodes[iW]->back().node ) + --nbP; // 1st and last points coincide + } + _bndNodes.resize( nbP ); + + // fill boundary points + BRepMesh::Array1OfVertexOfDelaun bndVert( 1, 1 + nbP ); + BRepMesh_Vertex v( 0, 0, BRepMesh_Frontier ); + for ( size_t iW = 0; iW < boundaryNodes.size(); ++iW ) + { + const UVPtStructVec& bndPnt = *boundaryNodes[iW]; + int i = 0, nb = bndPnt.size(); + if ( bndPnt[0].node == bndPnt.back().node ) + --nb; + for ( ; i < nb; ++i, ++iP ) + { + _bndNodes[ iP-1 ] = bndPnt[i].node; + bndPnt[i].node->setIsMarked( true ); + + v.ChangeCoord() = bndPnt[i].UV().Multiplied( _scale ); + bndVert( iP ) = v; + } + } + + // triangulate the srcFace in 2D + BRepMesh_Delaun Delaunay( bndVert ); + _triaDS = Delaunay.Result(); +} + +//================================================================================ +/*! + * \brief Prepare to the exploration of nodes + */ +//================================================================================ + +void SMESH_Delaunay::InitTraversal(const int nbNodesToVisit) +{ + _nbNodesToVisit = (size_t) nbNodesToVisit; + _nbVisitedNodes = _iBndNode = 0; + _noTriQueue.clear(); +} + +//================================================================================ +/*! + * \brief Return a node with its Barycentric Coordinates within the triangle + * defined by its node indices (zero based) + * \param [out] bc - Barycentric Coordinates of the returned node + * \param [out] triaNodes - indices of triangle nodes + * \return const SMDS_MeshNode* - the next node or NULL + */ +//================================================================================ + +const SMDS_MeshNode* SMESH_Delaunay::NextNode( double bc[3], int triaNodes[3] ) +{ + while ( _nbVisitedNodes < _nbNodesToVisit ) + { + while ( !_noTriQueue.empty() ) + { + const SMDS_MeshNode* node = _noTriQueue.front().first; + const BRepMesh_Triangle* tria = _noTriQueue.front().second; + _noTriQueue.pop_front(); + if ( node->isMarked() ) + continue; + ++_nbVisitedNodes; + node->setIsMarked( true ); + + // find a Delaunay triangle containing the src node + gp_XY uv = getNodeUV( _face, node ); + tria = FindTriangle( uv, tria, bc, triaNodes ); + if ( tria ) + { + addCloseNodes( node, tria, _faceID, _noTriQueue ); + return node; + } + } + for ( ; _iBndNode < _bndNodes.size() && _noTriQueue.empty(); ++_iBndNode ) + { + if ( const BRepMesh_Triangle* tria = GetTriangleNear( _iBndNode )) + addCloseNodes( _bndNodes[ _iBndNode ], tria, _faceID, _noTriQueue ); + } + if ( _noTriQueue.empty() ) + break; + } + + // if ( _nbVisitedNodes < _nbNodesToVisit ) + // _nbVisitedNodes = std::numeric_limits::max(); + return NULL; +} + +//================================================================================ +/*! + * \brief Find a Delaunay triangle containing a given 2D point and return + * barycentric coordinates within the found triangle + */ +//================================================================================ + +const BRepMesh_Triangle* SMESH_Delaunay::FindTriangle( const gp_XY& UV, + const BRepMesh_Triangle* tria, + double bc[3], + int triaNodes[3] ) +{ + int nodeIDs[3]; + gp_XY nodeUVs[3]; + int linkIDs[3]; + Standard_Boolean ori[3]; + + gp_XY uv = UV.Multiplied( _scale ); + + while ( tria ) + { + // check if the uv is in tria + + _triaDS->ElementNodes( *tria, nodeIDs ); + nodeUVs[0] = _triaDS->GetNode( nodeIDs[0] ).Coord(); + nodeUVs[1] = _triaDS->GetNode( nodeIDs[1] ).Coord(); + nodeUVs[2] = _triaDS->GetNode( nodeIDs[2] ).Coord(); + + if ( _triaDS->GetNode( nodeIDs[0] ).Movability() == BRepMesh_Frontier && + _triaDS->GetNode( nodeIDs[1] ).Movability() == BRepMesh_Frontier && + _triaDS->GetNode( nodeIDs[2] ).Movability() == BRepMesh_Frontier ) + { + SMESH_MeshAlgos::GetBarycentricCoords( uv, + nodeUVs[0], nodeUVs[1], nodeUVs[2], + bc[0], bc[1] ); + if ( bc[0] >= 0 && bc[1] >= 0 && bc[0] + bc[1] <= 1 ) + { + bc[2] = 1 - bc[0] - bc[1]; + triaNodes[0] = nodeIDs[0] - 1; + triaNodes[1] = nodeIDs[1] - 1; + triaNodes[2] = nodeIDs[2] - 1; + return tria; + } + } + + // look for a neighbor triangle, which is adjacent to a link intersected + // by a segment( triangle center -> uv ) + + gp_XY gc = ( nodeUVs[0] + nodeUVs[1] + nodeUVs[2] ) / 3.; + gp_XY seg = uv - gc; + + tria->Edges( linkIDs, ori ); + int triaID = _triaDS->IndexOf( *tria ); + tria = 0; + + for ( int i = 0; i < 3; ++i ) + { + const BRepMesh_PairOfIndex & triIDs = _triaDS->ElementsConnectedTo( linkIDs[i] ); + if ( triIDs.Extent() < 2 ) + continue; // no neighbor triangle + + // check if a link intersects gc2uv + const BRepMesh_Edge & link = _triaDS->GetLink( linkIDs[i] ); + const BRepMesh_Vertex & n1 = _triaDS->GetNode( link.FirstNode() ); + const BRepMesh_Vertex & n2 = _triaDS->GetNode( link.LastNode() ); + gp_XY uv1 = n1.Coord(); + gp_XY lin = n2.Coord() - uv1; // link direction + + double crossSegLin = seg ^ lin; + if ( Abs( crossSegLin ) < std::numeric_limits::min() ) + continue; // parallel + + double uSeg = ( uv1 - gc ) ^ lin / crossSegLin; + if ( 0. <= uSeg && uSeg <= 1. ) + { + tria = & _triaDS->GetElement( triIDs.Index( 1 + ( triIDs.Index(1) == triaID ))); + if ( tria->Movability() != BRepMesh_Deleted ) + break; + } + } + } + return tria; +} + +//================================================================================ +/*! + * \brief Return a triangle sharing a given boundary node + * \param [in] iBndNode - index of the boundary node + * \return const BRepMesh_Triangle* - a found triangle + */ +//================================================================================ + +const BRepMesh_Triangle* SMESH_Delaunay::GetTriangleNear( int iBndNode ) +{ + int nodeIDs[3]; + int nbNbNodes = _bndNodes.size(); + const BRepMesh::ListOfInteger & linkIds = _triaDS->LinksConnectedTo( iBndNode + 1 ); + BRepMesh::ListOfInteger::const_iterator iLink = linkIds.cbegin(); + for ( ; iLink != linkIds.cend(); ++iLink ) + { + const BRepMesh_PairOfIndex & triaIds = _triaDS->ElementsConnectedTo( *iLink ); + { + const BRepMesh_Triangle& tria = _triaDS->GetElement( triaIds.Index(1) ); + if ( tria.Movability() != BRepMesh_Deleted ) + { + _triaDS->ElementNodes( tria, nodeIDs ); + if ( nodeIDs[0]-1 < nbNbNodes && + nodeIDs[1]-1 < nbNbNodes && + nodeIDs[2]-1 < nbNbNodes ) + return &tria; + } + } + if ( triaIds.Extent() > 1 ) + { + const BRepMesh_Triangle& tria = _triaDS->GetElement( triaIds.Index(2) ); + if ( tria.Movability() != BRepMesh_Deleted ) + { + _triaDS->ElementNodes( tria, nodeIDs ); + if ( nodeIDs[0]-1 < nbNbNodes && + nodeIDs[1]-1 < nbNbNodes && + nodeIDs[2]-1 < nbNbNodes ) + return &tria; + } + } + } + return 0; +} + +//================================================================================ +/*! + * \brief Return UV of the i-th source boundary node (zero based) + */ +//================================================================================ + +gp_XY SMESH_Delaunay::GetBndUV(const int iNode) const +{ + return _triaDS->GetNode( iNode+1 ).Coord(); +} + +//================================================================================ +/*! + * \brief Add non-marked nodes surrounding a given one to a queue + */ +//================================================================================ + +void SMESH_Delaunay::addCloseNodes( const SMDS_MeshNode* node, + const BRepMesh_Triangle* tria, + const int faceID, + TNodeTriaList & _noTriQueue ) +{ + // find in-FACE nodes + SMDS_ElemIteratorPtr elems = node->GetInverseElementIterator(SMDSAbs_Face); + while ( elems->more() ) + { + const SMDS_MeshElement* elem = elems->next(); + if ( elem->getshapeId() == faceID ) + { + for ( int i = 0, nb = elem->NbNodes(); i < nb; ++i ) + { + const SMDS_MeshNode* n = elem->GetNode( i ); + if ( !n->isMarked() /*&& n->getshapeId() == faceID*/ ) + _noTriQueue.push_back( std::make_pair( n, tria )); + } + } + } +} + +//================================================================================ +/*! + * \brief Write a python script that creates an equal mesh in Mesh module + */ +//================================================================================ + +void SMESH_Delaunay::ToPython() const +{ + SMESH_Comment text; + text << "import salome, SMESH\n"; + text << "salome.salome_init()\n"; + text << "from salome.smesh import smeshBuilder\n"; + text << "smesh = smeshBuilder.New(salome.myStudy)\n"; + text << "mesh=smesh.Mesh()\n"; + const char* endl = "\n"; + + for ( int i = 0; i < _triaDS->NbNodes(); ++i ) + { + const BRepMesh_Vertex& v = _triaDS->GetNode( i+1 ); + text << "mesh.AddNode( " << v.Coord().X() << ", " << v.Coord().Y() << ", 0 )" << endl; + } + + int nodeIDs[3]; + for ( int i = 0; i < _triaDS->NbElements(); ++i ) + { + const BRepMesh_Triangle& t = _triaDS->GetElement( i+1 ); + if ( t.Movability() == BRepMesh_Deleted ) + continue; + _triaDS->ElementNodes( t, nodeIDs ); + text << "mesh.AddFace([ " << nodeIDs[0] << ", " << nodeIDs[1] << ", " << nodeIDs[2] << " ])" << endl; + } + + const char* fileName = "/tmp/Delaunay.py"; + SMESH_File file( fileName, false ); + file.remove(); + file.openForWriting(); + file.write( text.c_str(), text.size() ); + cout << "execfile( '" << fileName << "')" << endl; +} diff --git a/src/SMESHUtils/SMESH_Delaunay.hxx b/src/SMESHUtils/SMESH_Delaunay.hxx new file mode 100644 index 000000000..71707118b --- /dev/null +++ b/src/SMESHUtils/SMESH_Delaunay.hxx @@ -0,0 +1,115 @@ +// Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// File : SMESH_Delaunay.hxx +// Created : Wed Apr 19 15:42:54 2017 +// Author : Edward AGAPOV (eap) + + +#ifndef __SMESH_Delaunay_HXX__ +#define __SMESH_Delaunay_HXX__ + +#include "SMESH_TypeDefs.hxx" + +#include + +/*! + * \brief Create a Delaunay triangulation of nodes on a face boundary + * and provide exploration of nodes shared by elements lying on + * the face. For a returned node, also return a Delaunay triangle + * the node lies in and its Barycentric Coordinates within the triangle. + * Only non-marked nodes are visited. Boundary nodes given at the construction + * are not returned. + * + * For usage, this class needs to be subclassed to implement getNodeUV(); + */ +class SMESHUtils_EXPORT SMESH_Delaunay +{ + public: + + // construct a Delaunay triangulation of given boundary nodes + SMESH_Delaunay(const std::vector< const UVPtStructVec* > & boundaryNodes, + const TopoDS_Face& face, + const int faceID); + + virtual ~SMESH_Delaunay() {} + + // prepare to the exploration of nodes + void InitTraversal(const int nbNodesToVisit = -1); + + // return a node with its Barycentric Coordinates within the triangle + // defined by its node indices (zero based) + const SMDS_MeshNode* NextNode( double bc[3], int triaNodes[3] ); + + // return nb of nodes returned by NextNode() + int NbVisitedNodes() const { return _nbVisitedNodes; } + + + // find a triangle containing an UV, starting from a given triangle; + // return barycentric coordinates of the UV and the found triangle (indices are zero based). + const BRepMesh_Triangle* FindTriangle( const gp_XY& uv, + const BRepMesh_Triangle* bmTria, + double bc[3], + int triaNodes[3]); + + // return any Delaunay triangle neighboring a given boundary node (zero based) + const BRepMesh_Triangle* GetTriangleNear( int iBndNode ); + + // return source boundary nodes + const std::vector< const SMDS_MeshNode* >& GetBndNodes() const { return _bndNodes; } + + // return UV of the i-th source boundary node (zero based) + gp_XY GetBndUV(const int iNode) const; + + // return scale factor to convert real UV to/from UV used for Delauney meshing: + // delauney_UV = real_UV * scale + const gp_XY& GetScale() const { return _scale; } + + void ToPython() const; + + Handle(BRepMesh_DataStructureOfDelaun) GetDS() { return _triaDS; } + + protected: + + // container of a node and a triangle serving as a start while searching a + // triangle including the node UV + typedef std::list< std::pair< const SMDS_MeshNode*, const BRepMesh_Triangle* > > TNodeTriaList; + + // return UV of a node on the face + virtual gp_XY getNodeUV( const TopoDS_Face& face, const SMDS_MeshNode* node ) const = 0; + + // add non-marked nodes surrounding a given one to a queue + static void addCloseNodes( const SMDS_MeshNode* node, + const BRepMesh_Triangle* bmTria, + const int faceID, + TNodeTriaList & noTriQueue ); + + const TopoDS_Face& _face; + int _faceID; + std::vector< const SMDS_MeshNode* > _bndNodes; + gp_XY _scale; + Handle(BRepMesh_DataStructureOfDelaun) _triaDS; + size_t _nbNodesToVisit, _nbVisitedNodes, _iBndNode; + TNodeTriaList _noTriQueue; + +}; + +#endif diff --git a/src/SMESH_I/SMESH_Filter_i.cxx b/src/SMESH_I/SMESH_Filter_i.cxx index 739a4d904..9432b9117 100644 --- a/src/SMESH_I/SMESH_Filter_i.cxx +++ b/src/SMESH_I/SMESH_Filter_i.cxx @@ -4109,7 +4109,7 @@ static const char** getFunctNames() #ifdef _DEBUG_ // check if functName is complete, compilation failure means that enum FunctorType changed const int nbFunctors = sizeof(functName) / sizeof(const char*); - int _assert[( nbFunctors == SMESH::FT_Undefined + 1 ) ? 2 : -1 ]; _assert[0]=_assert[1]; + int _assert[( nbFunctors == SMESH::FT_Undefined + 1 ) ? 2 : -1 ]; _assert[0]=_assert[1]=0; #endif return functName; diff --git a/src/SMESH_I/SMESH_Gen_i.cxx b/src/SMESH_I/SMESH_Gen_i.cxx index 44cbff77a..916a1f222 100644 --- a/src/SMESH_I/SMESH_Gen_i.cxx +++ b/src/SMESH_I/SMESH_Gen_i.cxx @@ -2560,7 +2560,7 @@ SMESH_Gen_i::ConcatenateCommon(const SMESH::ListOfIDSources& theMeshesArray, const char* typeNames[] = { "All","Nodes","Edges","Faces","Volumes","0DElems","Balls" }; { // check of typeNames: compilation failure mains that NB_ELEMENT_TYPES changed: const int nbNames = sizeof(typeNames) / sizeof(const char*); - int _assert[( nbNames == SMESH::NB_ELEMENT_TYPES ) ? 2 : -1 ]; _assert[0]=_assert[1]; + int _assert[( nbNames == SMESH::NB_ELEMENT_TYPES ) ? 2 : -1 ]; _assert[0]=_assert[1]=0; } string groupName = "Gr"; SALOMEDS::SObject_wrap aMeshSObj = ObjectToSObject( myCurrentStudy, theMeshesArray[i] ); diff --git a/src/SMESH_I/SMESH_Mesh_i.cxx b/src/SMESH_I/SMESH_Mesh_i.cxx index a075dadce..268dfd829 100644 --- a/src/SMESH_I/SMESH_Mesh_i.cxx +++ b/src/SMESH_I/SMESH_Mesh_i.cxx @@ -4771,6 +4771,35 @@ CORBA::Long SMESH_Mesh_i::FindElementByNodes(const SMESH::long_array& nodes) return elemID; } +//================================================================================ +/*! + * \brief Return elements including all given nodes. + */ +//================================================================================ + +SMESH::long_array* SMESH_Mesh_i::GetElementsByNodes(const SMESH::long_array& nodes, + SMESH::ElementType elemType) +{ + if ( _preMeshInfo ) + _preMeshInfo->FullLoadFromFile(); + + SMESH::long_array_var result = new SMESH::long_array(); + + if ( SMESHDS_Mesh* mesh = _impl->GetMeshDS() ) + { + vector< const SMDS_MeshNode * > nn( nodes.length() ); + for ( CORBA::ULong i = 0; i < nodes.length(); ++i ) + nn[i] = mesh->FindNode( nodes[i] ); + + std::vector elems; + mesh->GetElementsByNodes( nn, elems, (SMDSAbs_ElementType) elemType ); + result->length( elems.size() ); + for ( size_t i = 0; i < elems.size(); ++i ) + result[i] = elems[i]->GetID(); + } + return result._retn(); +} + //============================================================================= /*! * Returns true if given element is polygon diff --git a/src/SMESH_I/SMESH_Mesh_i.hxx b/src/SMESH_I/SMESH_Mesh_i.hxx index 8aaa79bae..62d08763d 100644 --- a/src/SMESH_I/SMESH_Mesh_i.hxx +++ b/src/SMESH_I/SMESH_Mesh_i.hxx @@ -536,8 +536,8 @@ public: * Returns true if given node is medium node * in one of quadratic elements */ - CORBA::Boolean IsMediumNodeOfAnyElem(CORBA::Long idn, - SMESH::ElementType theElemType); + CORBA::Boolean IsMediumNodeOfAnyElem(CORBA::Long idn, + SMESH::ElementType elemType); /*! * Returns number of edges for given element @@ -563,6 +563,12 @@ public: */ CORBA::Long FindElementByNodes(const SMESH::long_array& nodes); + /*! + * Return elements including all given nodes. + */ + SMESH::long_array* GetElementsByNodes(const SMESH::long_array& nodes, + SMESH::ElementType elemType); + /*! * Returns true if given element is polygon */ diff --git a/src/SMESH_SWIG/CMakeLists.txt b/src/SMESH_SWIG/CMakeLists.txt index 057506db6..b57d2b4df 100644 --- a/src/SMESH_SWIG/CMakeLists.txt +++ b/src/SMESH_SWIG/CMakeLists.txt @@ -22,8 +22,6 @@ # scripts / static SET(_bin_SCRIPTS smesh.py - batchmode_smesh.py - batchmode_mefisto.py ex00_all.py ex01_cube2build.py ex02_cube2primitive.py diff --git a/src/SMESH_SWIG/batchmode_mefisto.py b/src/SMESH_SWIG/batchmode_mefisto.py deleted file mode 100644 index baf57efa9..000000000 --- a/src/SMESH_SWIG/batchmode_mefisto.py +++ /dev/null @@ -1,127 +0,0 @@ -# -*- coding: iso-8859-1 -*- -# Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE -# -# Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -# CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -# -# This library is free software; you can redistribute it and/or -# modify it under the terms of the GNU Lesser General Public -# License as published by the Free Software Foundation; either -# version 2.1 of the License, or (at your option) any later version. -# -# This library is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# Lesser General Public License for more details. -# -# You should have received a copy of the GNU Lesser General Public -# License along with this library; if not, write to the Free Software -# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA -# -# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com -# - -import os -import re - -import batchmode_salome -import batchmode_geompy -import batchmode_smesh -from salome.StdMeshers import StdMeshersBuilder - -smesh = batchmode_smesh.smesh -smesh.SetCurrentStudy(batchmode_salome.myStudy) - -def CreateMesh (theFileName, area, len = None, nbseg = None): - - if not(os.path.isfile(theFileName)) or re.search("\.brep$", theFileName) is None : - print "Incorrect file name !" - return - - if (len is None) and (nbseg is None): - print "Define length or number of segments !" - return - - if (len is not None) and (nbseg is not None): - print "Only one Hypothesis (from length and number of segments) can be defined !" - return - - - # ---- Import shape from BREP file and add it to the study - shape_mesh = batchmode_geompy.Import(theFileName, "BREP") - Id_shape = batchmode_geompy.addToStudy(shape_mesh, "shape_mesh") - - - # ---- SMESH - print "-------------------------- create mesh" - mesh = smesh.Mesh(shape_mesh) - - print "-------------------------- create Hypothesis" - if (len is not None): - print "-------------------------- LocalLength" - algoReg = mesh.Segment() - hypLength1 = algoReg.LocalLength(len) - print "Hypothesis type : ", hypLength1.GetName() - print "Hypothesis ID : ", hypLength1.GetId() - print "Hypothesis Value: ", hypLength1.GetLength() - - if (nbseg is not None): - print "-------------------------- NumberOfSegments" - algoReg = mesh.Segment() - hypNbSeg1 = algoReg.NumberOfSegments(nbseg) - print "Hypothesis type : ", hypNbSeg1.GetName() - print "Hypothesis ID : ", hypNbSeg1.GetId() - print "Hypothesis Value: ", hypNbSeg1.GetNumberOfSegments() - - if (area == "LengthFromEdges"): - print "-------------------------- LengthFromEdges" - algoMef = mesh.Triangle() - hypLengthFromEdges = algoMef.LengthFromEdges(1) - print "Hypothesis type : ", hypLengthFromEdges.GetName() - print "Hypothesis ID : ", hypLengthFromEdges.GetId() - print "LengthFromEdges Mode: ", hypLengthFromEdges.GetMode() - - else: - print "-------------------------- MaxElementArea" - algoMef = mesh.Triangle() - hypArea1 = algoMef.MaxElementArea(area) - print "Hypothesis type : ", hypArea1.GetName() - print "Hypothesis ID : ", hypArea1.GetId() - print "Hypothesis Value: ", hypArea1.GetMaxElementArea() - - - print "-------------------------- Regular_1D" - listHyp = algoReg.GetCompatibleHypothesis() - for hyp in listHyp: - print hyp - - print "Algo name: ", algoReg.GetName() - print "Algo ID : ", algoReg.GetId() - - print "-------------------------- MEFISTO_2D" - listHyp = algoMef.GetCompatibleHypothesis() - for hyp in listHyp: - print hyp - - print "Algo name: ", algoMef.GetName() - print "Algo ID : ", algoMef.GetId() - - - # ---- add hypothesis to shape - - print "-------------------------- compute mesh" - ret = mesh.Compute() - print "Compute Mesh .... ", - print ret - log = mesh.GetLog(0); # no erase trace - #for linelog in log: - # print linelog - - print "------------ INFORMATION ABOUT MESH ------------" - - print "Number of nodes : ", mesh.NbNodes() - print "Number of edges : ", mesh.NbEdges() - print "Number of faces : ", mesh.NbFaces() - print "Number of triangles: ", mesh.NbTriangles() - - return mesh diff --git a/src/SMESH_SWIG/batchmode_smesh.py b/src/SMESH_SWIG/batchmode_smesh.py deleted file mode 100644 index c31f82a04..000000000 --- a/src/SMESH_SWIG/batchmode_smesh.py +++ /dev/null @@ -1,324 +0,0 @@ -# -*- coding: iso-8859-1 -*- -# Copyright (C) 2007-2016 CEA/DEN, EDF R&D, OPEN CASCADE -# -# Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -# CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -# -# This library is free software; you can redistribute it and/or -# modify it under the terms of the GNU Lesser General Public -# License as published by the Free Software Foundation; either -# version 2.1 of the License, or (at your option) any later version. -# -# This library is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# Lesser General Public License for more details. -# -# You should have received a copy of the GNU Lesser General Public -# License along with this library; if not, write to the Free Software -# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA -# -# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com -# - -# File : batchmode_smesh.py -# Author : Oksana TCHEBANOVA -# Module : SMESH -# $Header$ -# -from batchmode_salome import * -from batchmode_geompy import ShapeType -import SMESH - -#-------------------------------------------------------------------------- -modulecatalog = naming_service.Resolve("/Kernel/ModulCatalog") - -smesh = lcc.FindOrLoadComponent("FactoryServer", "SMESH") -smesh.SetCurrentStudy(myStudy) -myStudyBuilder = myStudy.NewBuilder() - -if myStudyBuilder is None: - raise RuntimeError, " Null myStudyBuilder" - -father = myStudy.FindComponent("SMESH") -if father is None: - father = myStudyBuilder.NewComponent("SMESH") - FName = myStudyBuilder.FindOrCreateAttribute(father, "AttributeName") - Comp = modulecatalog.GetComponent("SMESH") - FName.SetValue(Comp._get_componentusername()) - aPixmap = myStudyBuilder.FindOrCreateAttribute(father, "AttributePixMap") - aPixmap.SetPixMap("ICON_OBJBROWSER_Mesh") - -myStudyBuilder.DefineComponentInstance(father,smesh) - -mySComponentMesh = father._narrow(SALOMEDS.SComponent) - -Tag_HypothesisRoot = 1 -Tag_AlgorithmsRoot = 2 - -Tag_RefOnShape = 1 -Tag_RefOnAppliedHypothesis = 2 -Tag_RefOnAppliedAlgorithms = 3 - -Tag_SubMeshOnVertex = 4 -Tag_SubMeshOnEdge = 5 -Tag_SubMeshOnFace = 6 -Tag_SubMeshOnSolid = 7 -Tag_SubMeshOnCompound = 8 - -Tag = {"HypothesisRoot":1,"AlgorithmsRoot":2,"RefOnShape":1,"RefOnAppliedHypothesis":2,"RefOnAppliedAlgorithms":3,"SubMeshOnVertex":4,"SubMeshOnEdge":5,"SubMeshOnFace":6,"SubMeshOnSolid":7,"SubMeshOnCompound":8} - -#------------------------------------------------------------ -def Init(): - pass -#------------------------------------------------------------ -def AddNewMesh(IOR): - # VSR: added temporarily - objects are published automatically by the engine - aSO = myStudy.FindObjectIOR( IOR ) - if aSO is not None: - return aSO.GetID() - # VSR ###################################################################### - - res,HypothesisRoot = mySComponentMesh.FindSubObject ( Tag_HypothesisRoot ) - if HypothesisRoot is None or res == 0: - HypothesisRoot = myStudyBuilder.NewObjectToTag(mySComponentMesh, Tag_HypothesisRoot) - aName = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributeName") - aName.SetValue("Hypotheses") - aPixmap = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_hypo.png" ) - aSelAttr = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributeSelectable") - aSelAttr.SetSelectable(0) - - res, AlgorithmsRoot = mySComponentMesh.FindSubObject (Tag_AlgorithmsRoot) - if AlgorithmsRoot is None or res == 0: - AlgorithmsRoot = myStudyBuilder.NewObjectToTag (mySComponentMesh, Tag_AlgorithmsRoot) - aName = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributeName") - aName.SetValue("Algorithms") - aPixmap = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_algo.png" ) - aSelAttr = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributeSelectable") - aSelAttr.SetSelectable(0) - - HypothesisRoot = HypothesisRoot._narrow(SALOMEDS.SObject) - newMesh = myStudyBuilder.NewObject(mySComponentMesh) - aPixmap = myStudyBuilder.FindOrCreateAttribute(newMesh, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_mesh.png" ) - anIOR = myStudyBuilder.FindOrCreateAttribute(newMesh, "AttributeIOR") - anIOR.SetValue(IOR) - return newMesh.GetID() - -#------------------------------------------------------------ -def AddNewHypothesis(IOR): - # VSR: added temporarily - objects are published automatically by the engine - aSO = myStudy.FindObjectIOR( IOR ) - if aSO is not None: - return aSO.GetID() - # VSR ###################################################################### - - res, HypothesisRoot = mySComponentMesh.FindSubObject (Tag_HypothesisRoot) - if HypothesisRoot is None or res == 0: - HypothesisRoot = myStudyBuilder.NewObjectToTag (mySComponentMesh, Tag_HypothesisRoot) - aName = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributeName") - aName.SetValue("Hypotheses") - aSelAttr = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributeSelectable") - aSelAttr.SetSelectable(0) - aPixmap = myStudyBuilder.FindOrCreateAttribute(HypothesisRoot, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_hypo.png" ) - - # Add New Hypothesis - newHypo = myStudyBuilder.NewObject(HypothesisRoot) - aPixmap = myStudyBuilder.FindOrCreateAttribute(newHypo, "AttributePixMap") - H = orb.string_to_object(IOR) - aType = H.GetName() - aPixmap.SetPixMap( "mesh_tree_hypo.png_" + aType ) - anIOR = myStudyBuilder.FindOrCreateAttribute(newHypo, "AttributeIOR") - anIOR.SetValue(IOR) - return newHypo.GetID() - -#------------------------------------------------------------ -def AddNewAlgorithms(IOR): - # VSR: added temporarily - objects are published automatically by the engine - aSO = myStudy.FindObjectIOR( IOR ) - if aSO is not None: - return aSO.GetID() - # VSR ###################################################################### - - res, AlgorithmsRoot = mySComponentMesh.FindSubObject (Tag_AlgorithmsRoot) - if AlgorithmsRoot is None or res == 0: - AlgorithmsRoot = myStudyBuilde.NewObjectToTag (mySComponentMesh, Tag_AlgorithmsRoot) - aName = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributeName") - aName.SetValue("Algorithms") - aSelAttr = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributeSelectable") - aSelAttr.SetSelectable(0) - aPixmap = myStudyBuilder.FindOrCreateAttribute(AlgorithmsRoot, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_algo.png" ) - - # Add New Algorithms - newHypo = myStudyBuilder.NewObject(AlgorithmsRoot) - aPixmap = myStudyBuilder.FindOrCreateAttribute(newHypo, "AttributePixMap") - aPixmap = anAttr._narrow(SALOMEDS.AttributePixMap) - H = orb.string_to_object(IOR) - aType = H.GetName(); #QString in fact - aPixmap.SetPixMap( "mesh_tree_algo.png_" + aType ) - anIOR = myStudyBuilder.FindOrCreateAttribute(newHypo, "AttributeIOR") - anIOR.SetValue(IOR) - return newHypo.GetID() - - -#------------------------------------------------------------ -def SetShape(ShapeEntry, MeshEntry): - SO_MorSM = myStudy.FindObjectID( MeshEntry ) - SO_GeomShape = myStudy.FindObjectID( ShapeEntry ) - - if SO_MorSM is not None and SO_GeomShape is not None : - # VSR: added temporarily - shape reference is published automatically by the engine - res, Ref = SO_MorSM.FindSubObject( Tag_RefOnShape ) - if res == 1 : - return - # VSR ###################################################################### - - SO = myStudyBuilder.NewObjectToTag (SO_MorSM, Tag_RefOnShape) - myStudyBuilder.Addreference (SO,SO_GeomShape) - - -#------------------------------------------------------------ -def SetHypothesis(Mesh_Or_SubMesh_Entry, Hypothesis_Entry): - SO_MorSM = myStudy.FindObjectID( Mesh_Or_SubMesh_Entry ) - SO_Hypothesis = myStudy.FindObjectID( Hypothesis_Entry ) - - if SO_MorSM is not None and SO_Hypothesis is not None : - - #Find or Create Applied Hypothesis root - res, AHR = SO_MorSM.FindSubObject (Tag_RefOnAppliedHypothesis) - if AHR is None or res == 0: - AHR = myStudyBuilder.NewObjectToTag (SO_MorSM, Tag_RefOnAppliedHypothesis) - aName = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributeName") - - # The same name as in SMESH_Mesh_i::AddHypothesis() ################## - aName.SetValue("Applied hypotheses") - - aSelAttr = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributeSelectable") - aSelAttr.SetSelectable(0) - aPixmap = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_hypo.png" ) - - # VSR: added temporarily - reference to applied hypothesis is published automatically by the engine - else : - it = myStudy.NewChildIterator(AHR) - while it.More() : - res, Ref = it.Value().ReferencedObject() - if res and Ref is not None and Ref.GetID() == Hypothesis_Entry : - return - it.Next() - # VSR ###################################################################### - - SO = myStudyBuilder.NewObject(AHR) - myStudyBuilder.Addreference (SO,SO_Hypothesis) - -#------------------------------------------------------------ -def SetAlgorithms(Mesh_Or_SubMesh_Entry, Algorithms_Entry): - SO_MorSM = myStudy.FindObjectID( Mesh_Or_SubMesh_Entry ) - SO_Algorithms = myStudy.FindObjectID( Algorithms_Entry ) - if SO_MorSM != None and SO_Algorithms != None : - #Find or Create Applied Algorithms root - res, AHR = SO_MorSM.FindSubObject (Tag_RefOnAppliedAlgorithms) - if AHR is None or res == 0: - AHR = myStudyBuilder.NewObjectToTag (SO_MorSM, Tag_RefOnAppliedAlgorithms) - aName = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributeName") - - # The same name as in SMESH_Mesh_i::AddHypothesis() ################## - aName.SetValue("Applied algorithms") - - aSelAttr = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributeSelectable") - aSelAttr.SetSelectable(0) - aPixmap = myStudyBuilder.FindOrCreateAttribute(AHR, "AttributePixMap") - aPixmap.SetPixMap( "mesh_tree_algo.png" ) - - # VSR: added temporarily - reference to applied hypothesis is published automatically by the engine - else : - it = myStudy.NewChildIterator(AHR) - while it.More() : - res, Ref = it.Value().ReferencedObject() - if res and Ref is not None and Ref.GetID() == Algorithms_Entry : - return - it.Next() - # VSR ###################################################################### - - SO = myStudyBuilder.NewObject(AHR) - myStudyBuilder.Addreference (SO,SO_Algorithms) - - -#------------------------------------------------------------ -def UnSetHypothesis( Applied_Hypothesis_Entry ): - SO_Applied_Hypothesis = myStudy.FindObjectID( Applied_Hypothesis_Entry ) - if SO_Applied_Hypothesis : - myStudyBuilder.RemoveObject(SO_Applied_Hypothesis) - - -#------------------------------------------------------------ -def AddSubMesh ( SO_Mesh_Entry, SM_IOR, ST): - # VSR: added temporarily - objects are published automatically by the engine - aSO = myStudy.FindObjectIOR( SM_IOR ) - if aSO is not None: - return aSO.GetID() - # VSR ###################################################################### - - SO_Mesh = myStudy.FindObjectID( SO_Mesh_Entry ) - if ( SO_Mesh ) : - - if ST == ShapeType["COMPSOLID"] : - Tag_Shape = Tag_SubMeshOnSolid - Name = "SubMeshes on Solid" - elif ST == ShapeType["FACE"] : - Tag_Shape = Tag_SubMeshOnFace - Name = "SubMeshes on Face" - elif ST == ShapeType["EDGE"] : - Tag_Shape = Tag_SubMeshOnEdge - Name = "SubMeshes on Edge" - elif ST == ShapeType["VERTEX"] : - Tag_Shape = Tag_SubMeshOnVertex - Name = "SubMeshes on Vertex" - else : - Tag_Shape = Tag_SubMeshOnCompound - Name = "SubMeshes on Compound" - - res, SubmeshesRoot = SO_Mesh.FindSubObject (Tag_Shape) - if SubmeshesRoot is None or res == 0: - SubmeshesRoot = myStudyBuilder.NewObjectToTag (SO_Mesh, Tag_Shape) - aName = myStudyBuilder.FindOrCreateAttribute(SubmeshesRoot, "AttributeName") - aName.SetValue(Name) - aSelAttr = myStudyBuilder.FindOrCreateAttribute(SubmeshesRoot, "AttributeSelectable") - aSelAttr.SetSelectable(0) - - SO = myStudyBuilder.NewObject (SubmeshesRoot) - anIOR = myStudyBuilder.FindOrCreateAttribute(SO, "AttributeIOR") - anIOR.SetValue(SM_IOR) - return SO.GetID() - - return None - -#------------------------------------------------------------ -def AddSubMeshOnShape (Mesh_Entry, GeomShape_Entry, SM_IOR, ST) : - # VSR: added temporarily - objects are published automatically by the engine - aSO = myStudy.FindObjectIOR( SM_IOR ) - if aSO is not None: - return aSO.GetID() - # VSR ###################################################################### - SO_GeomShape = myStudy.FindObjectID( GeomShape_Entry ) - if SO_GeomShape != None : - SM_Entry = AddSubMesh (Mesh_Entry,SM_IOR,ST) - SO_SM = myStudy.FindObjectID( SM_Entry ) - - if SO_SM != None : - SetShape (GeomShape_Entry, SM_Entry) - return SM_Entry - - return None - - -#------------------------------------------------------------ -def SetName(Entry, Name): - SO = myStudy.FindObjectID( Entry ) - if SO != None : - aName = myStudyBuilder.FindOrCreateAttribute(SO, "AttributeName") - aName.SetValue(Name) diff --git a/src/SMESH_SWIG/ex29_refine.py b/src/SMESH_SWIG/ex29_refine.py index 938d0eac1..72e0f631a 100644 --- a/src/SMESH_SWIG/ex29_refine.py +++ b/src/SMESH_SWIG/ex29_refine.py @@ -36,8 +36,11 @@ import os # Values # ------ +tmpDir = os.getenv('SALOME_TMP_DIR', '/tmp') +print "Output directory:", tmpDir + # Path for ".med" files -path = "/tmp/ex29_%s_" % os.getenv('USER','unknown') +path = os.path.join( tmpDir, "ex29_%s_" % os.getenv('USER','unknown')) # Name of the shape and the mesh name = "Carre" diff --git a/src/SMESH_SWIG/smeshBuilder.py b/src/SMESH_SWIG/smeshBuilder.py index d87d0d553..88d2b0b1a 100644 --- a/src/SMESH_SWIG/smeshBuilder.py +++ b/src/SMESH_SWIG/smeshBuilder.py @@ -2687,9 +2687,14 @@ class Mesh: ## Return an element based on all given nodes. # @ingroup l1_meshinfo - def FindElementByNodes(self,nodes): + def FindElementByNodes(self, nodes): return self.mesh.FindElementByNodes(nodes) + ## Return elements including all given nodes. + # @ingroup l1_meshinfo + def GetElementsByNodes(self, nodes, elemType=SMESH.ALL): + return self.mesh.GetElementsByNodes( nodes, elemType ) + ## Return true if the given element is a polygon # @ingroup l1_meshinfo def IsPoly(self, id): @@ -5177,10 +5182,11 @@ omniORB.registerObjref(SMESH._objref_SMESH_Pattern._NP_RepositoryId, Pattern) ## Private class used to bind methods creating algorithms to the class Mesh # class algoCreator: - def __init__(self): + def __init__(self, method): self.mesh = None self.defaultAlgoType = "" self.algoTypeToClass = {} + self.method = method # Store a python class of algorithm def add(self, algoClass): @@ -5194,25 +5200,48 @@ class algoCreator: # Create a copy of self and assign mesh to the copy def copy(self, mesh): - other = algoCreator() + other = algoCreator( self.method ) other.defaultAlgoType = self.defaultAlgoType - other.algoTypeToClass = self.algoTypeToClass + other.algoTypeToClass = self.algoTypeToClass other.mesh = mesh return other # Create an instance of algorithm def __call__(self,algo="",geom=0,*args): - algoType = self.defaultAlgoType - for arg in args + (algo,geom): - if isinstance( arg, geomBuilder.GEOM._objref_GEOM_Object ): - geom = arg - if isinstance( arg, str ) and arg: + algoType = "" + shape = 0 + if isinstance( algo, str ): + algoType = algo + elif ( isinstance( algo, geomBuilder.GEOM._objref_GEOM_Object ) and \ + not isinstance( geom, geomBuilder.GEOM._objref_GEOM_Object )): + shape = algo + elif algo: + args += (algo,) + + if isinstance( geom, geomBuilder.GEOM._objref_GEOM_Object ): + shape = geom + elif not algoType and isinstance( geom, str ): + algoType = geom + elif geom: + args += (geom,) + for arg in args: + if isinstance( arg, geomBuilder.GEOM._objref_GEOM_Object ) and not shape: + shape = arg + elif isinstance( arg, str ) and not algoType: algoType = arg + else: + import traceback, sys + msg = "Warning. Unexpected argument in mesh.%s() ---> %s" % ( self.method, arg ) + sys.stderr.write( msg + '\n' ) + tb = traceback.extract_stack(None,2) + traceback.print_list( [tb[0]] ) + if not algoType: + algoType = self.defaultAlgoType if not algoType and self.algoTypeToClass: algoType = self.algoTypeToClass.keys()[0] if self.algoTypeToClass.has_key( algoType ): #print "Create algo",algoType - return self.algoTypeToClass[ algoType ]( self.mesh, geom ) + return self.algoTypeToClass[ algoType ]( self.mesh, shape ) raise RuntimeError, "No class found for algo type %s" % algoType return None @@ -5294,7 +5323,7 @@ for pluginName in os.environ[ "SMESH_MeshersList" ].split( ":" ): if type( algo ).__name__ == 'classobj' and hasattr( algo, "meshMethod" ): #print " meshMethod:" , str(algo.meshMethod) if not hasattr( Mesh, algo.meshMethod ): - setattr( Mesh, algo.meshMethod, algoCreator() ) + setattr( Mesh, algo.meshMethod, algoCreator( algo.meshMethod )) pass getattr( Mesh, algo.meshMethod ).add( algo ) pass diff --git a/src/StdMeshers/StdMeshers_FaceSide.cxx b/src/StdMeshers/StdMeshers_FaceSide.cxx index 9f5929e36..b0d36a59c 100644 --- a/src/StdMeshers/StdMeshers_FaceSide.cxx +++ b/src/StdMeshers/StdMeshers_FaceSide.cxx @@ -88,13 +88,13 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, */ //================================================================================ -StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, - std::list& theEdges, - SMESH_Mesh* theMesh, - const bool theIsForward, - const bool theIgnoreMediumNodes, - SMESH_MesherHelper* theFaceHelper, - SMESH_ProxyMesh::Ptr theProxyMesh) +StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, + const std::list& theEdges, + SMESH_Mesh* theMesh, + const bool theIsForward, + const bool theIgnoreMediumNodes, + SMESH_MesherHelper* theFaceHelper, + SMESH_ProxyMesh::Ptr theProxyMesh) { int nbEdges = theEdges.size(); myEdge.resize ( nbEdges ); @@ -125,7 +125,7 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, SMESHDS_Mesh* meshDS = myProxyMesh->GetMeshDS(); int nbDegen = 0; - std::list::iterator edge = theEdges.begin(); + std::list::const_iterator edge = theEdges.begin(); for ( int index = 0; edge != theEdges.end(); ++index, ++edge ) { int i = theIsForward ? index : nbEdges-index-1; @@ -166,6 +166,14 @@ StdMeshers_FaceSide::StdMeshers_FaceSide(const TopoDS_Face& theFace, myC3dAdaptor[i].Load( C3d, 0, 0.5 * BRep_Tool::Tolerance( V )); } } + else if ( myEdgeLength[i] > DBL_MIN ) + { + Handle(Geom_Curve) C3d = BRep_Tool::Curve(myEdge[i],myFirst[i], myLast[i] ); + myC3dAdaptor[i].Load( C3d, myFirst[i], myLast[i] ); + if ( myEdge[i].Orientation() == TopAbs_REVERSED ) + std::swap( myFirst[i], myLast[i] ); + } + // reverse a proxy sub-mesh if ( !theIsForward ) reverseProxySubmesh( myEdge[i] ); @@ -1243,7 +1251,7 @@ gp_Pnt StdMeshers_FaceSide::Value3d(double U) const // check parametrization of curve if( !myIsUniform[i] ) { - double aLen3dU = r * myEdgeLength[i] * ( myFirst[i]>myLast[i] ? -1. : 1.); + double aLen3dU = r * myEdgeLength[i] * ( myFirst[i] > myLast[i] ? -1. : 1. ); GCPnts_AbscissaPoint AbPnt ( const_cast( myC3dAdaptor[i]), aLen3dU, myFirst[i] ); if( AbPnt.IsDone() ) { diff --git a/src/StdMeshers/StdMeshers_FaceSide.hxx b/src/StdMeshers/StdMeshers_FaceSide.hxx index 6759f50a9..1964638fe 100644 --- a/src/StdMeshers/StdMeshers_FaceSide.hxx +++ b/src/StdMeshers/StdMeshers_FaceSide.hxx @@ -82,13 +82,13 @@ public: /*! * \brief Wrap several edges. Edges must be properly ordered and oriented. */ - StdMeshers_FaceSide(const TopoDS_Face& theFace, - std::list& theEdges, - SMESH_Mesh* theMesh, - const bool theIsForward, - const bool theIgnoreMediumNodes, - SMESH_MesherHelper* theFaceHelper = NULL, - SMESH_ProxyMesh::Ptr theProxyMesh = SMESH_ProxyMesh::Ptr()); + StdMeshers_FaceSide(const TopoDS_Face& theFace, + const std::list& theEdges, + SMESH_Mesh* theMesh, + const bool theIsForward, + const bool theIgnoreMediumNodes, + SMESH_MesherHelper* theFaceHelper = NULL, + SMESH_ProxyMesh::Ptr theProxyMesh = SMESH_ProxyMesh::Ptr()); /*! * \brief Simulate a side from a vertex using data from other FaceSide */ @@ -120,13 +120,13 @@ public: { return StdMeshers_FaceSidePtr ( new StdMeshers_FaceSide( Face,Edge,Mesh,IsForward,IgnoreMediumNodes,FaceHelper,ProxyMesh )); } - static StdMeshers_FaceSidePtr New (const TopoDS_Face& Face, - std::list& Edges, - SMESH_Mesh* Mesh, - const bool IsForward, - const bool IgnoreMediumNodes, - SMESH_MesherHelper* FaceHelper = NULL, - SMESH_ProxyMesh::Ptr ProxyMesh = SMESH_ProxyMesh::Ptr()) + static StdMeshers_FaceSidePtr New (const TopoDS_Face& Face, + const std::list& Edges, + SMESH_Mesh* Mesh, + const bool IsForward, + const bool IgnoreMediumNodes, + SMESH_MesherHelper* FaceHelper = NULL, + SMESH_ProxyMesh::Ptr ProxyMesh = SMESH_ProxyMesh::Ptr()) { return StdMeshers_FaceSidePtr ( new StdMeshers_FaceSide( Face,Edges,Mesh,IsForward,IgnoreMediumNodes,FaceHelper,ProxyMesh )); } @@ -141,9 +141,11 @@ public: ( new StdMeshers_FaceSide( Side,Node,Pnt2d1,Pnt2d2,C2d,UFirst,ULast )); } static StdMeshers_FaceSidePtr New (UVPtStructVec& theSideNodes, - const TopoDS_Face& theFace = TopoDS_Face()) + const TopoDS_Face& theFace = TopoDS_Face(), + const TopoDS_Edge& theEdge = TopoDS_Edge(), + SMESH_Mesh* theMesh = 0) { - return StdMeshers_FaceSidePtr( new StdMeshers_FaceSide( theSideNodes, theFace )); + return StdMeshers_FaceSidePtr( new StdMeshers_FaceSide( theSideNodes, theFace, theEdge, theMesh )); } /*! diff --git a/src/StdMeshers/StdMeshers_Hexa_3D.hxx b/src/StdMeshers/StdMeshers_Hexa_3D.hxx index 3673cdba9..8b50468b7 100644 --- a/src/StdMeshers/StdMeshers_Hexa_3D.hxx +++ b/src/StdMeshers/StdMeshers_Hexa_3D.hxx @@ -54,6 +54,10 @@ public: virtual bool Evaluate(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape, MapShapeNbElems& aResMap); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); protected: diff --git a/src/StdMeshers/StdMeshers_Prism_3D.cxx b/src/StdMeshers/StdMeshers_Prism_3D.cxx index a7dae0283..ae8ccc430 100644 --- a/src/StdMeshers/StdMeshers_Prism_3D.cxx +++ b/src/StdMeshers/StdMeshers_Prism_3D.cxx @@ -51,6 +51,7 @@ #include #include #include +#include #include #include #include @@ -528,6 +529,27 @@ namespace { return nbSides; } + //================================================================================ + /*! + * \brief Set/get wire index to FaceQuadStruct + */ + //================================================================================ + + void setWireIndex( TFaceQuadStructPtr& quad, int iWire ) + { + quad->iSize = iWire; + } + int getWireIndex( const TFaceQuadStructPtr& quad ) + { + return quad->iSize; + } + + //================================================================================ + /*! + * \brief Print Python commands adding given points to a mesh + */ + //================================================================================ + void pointsToPython(const std::vector& p) { #ifdef _DEBUG_ @@ -558,6 +580,7 @@ StdMeshers_Prism_3D::StdMeshers_Prism_3D(int hypId, int studyId, SMESH_Gen* gen) //myProjectTriangles = false; mySetErrorToSM = true; // to pass an error to a sub-mesh of a current solid or not + myPrevBottomSM = 0; // last treated bottom sub-mesh with a suitable algorithm } //================================================================================ @@ -580,39 +603,6 @@ bool StdMeshers_Prism_3D::CheckHypothesis(SMESH_Mesh& a const TopoDS_Shape& aShape, SMESH_Hypothesis::Hypothesis_Status& aStatus) { - // Check shape geometry -/* PAL16229 - aStatus = SMESH_Hypothesis::HYP_BAD_GEOMETRY; - - // find not quadrangle faces - list< TopoDS_Shape > notQuadFaces; - int nbEdge, nbWire, nbFace = 0; - TopExp_Explorer exp( aShape, TopAbs_FACE ); - for ( ; exp.More(); exp.Next() ) { - ++nbFace; - const TopoDS_Shape& face = exp.Current(); - nbEdge = NSProjUtils::Count( face, TopAbs_EDGE, 0 ); - nbWire = NSProjUtils::Count( face, TopAbs_WIRE, 0 ); - if ( nbEdge!= 4 || nbWire!= 1 ) { - if ( !notQuadFaces.empty() ) { - if ( NSProjUtils::Count( notQuadFaces.back(), TopAbs_EDGE, 0 ) != nbEdge || - NSProjUtils::Count( notQuadFaces.back(), TopAbs_WIRE, 0 ) != nbWire ) - RETURN_BAD_RESULT("Different not quad faces"); - } - notQuadFaces.push_back( face ); - } - } - if ( !notQuadFaces.empty() ) - { - if ( notQuadFaces.size() != 2 ) - RETURN_BAD_RESULT("Bad nb not quad faces: " << notQuadFaces.size()); - - // check total nb faces - nbEdge = NSProjUtils::Count( notQuadFaces.back(), TopAbs_EDGE, 0 ); - if ( nbFace != nbEdge + 2 ) - RETURN_BAD_RESULT("Bad nb of faces: " << nbFace << " but must be " << nbEdge + 2); - } -*/ // no hypothesis aStatus = SMESH_Hypothesis::HYP_OK; return true; @@ -627,6 +617,7 @@ bool StdMeshers_Prism_3D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& theSh { SMESH_MesherHelper helper( theMesh ); myHelper = &helper; + myPrevBottomSM = 0; int nbSolids = helper.Count( theShape, TopAbs_SOLID, /*skipSame=*/false ); if ( nbSolids < 1 ) @@ -943,7 +934,7 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism, list< TopoDS_Edge >::iterator edge = thePrism.myBottomEdges.begin(); std::list< int >::iterator nbE = thePrism.myNbEdgesInWires.begin(); std::list< int > nbQuadsPerWire; - int iE = 0; + int iE = 0, iWire = 0; while ( edge != thePrism.myBottomEdges.end() ) { ++iE; @@ -977,7 +968,10 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism, return toSM( error(TCom("Composite 'horizontal' edges are not supported"))); } if ( faceMap.Add( face )) + { + setWireIndex( quadList.back(), iWire ); // for use in makeQuadsForOutInProjection() thePrism.myWallQuads.push_back( quadList ); + } break; } } @@ -995,6 +989,7 @@ bool StdMeshers_Prism_3D::getWallFaces( Prism_3D::TPrismTopo & thePrism, if ( iE == *nbE ) { iE = 0; + ++iWire; ++nbE; int nbQuadPrev = std::accumulate( nbQuadsPerWire.begin(), nbQuadsPerWire.end(), 0 ); nbQuadsPerWire.push_back( thePrism.myWallQuads.size() - nbQuadPrev ); @@ -1129,13 +1124,8 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism) return toSM( error( SMESH_ComputeError::New(COMPERR_CANCELED))); // Assure the bottom is meshed - SMESH_subMesh * botSM = myHelper->GetMesh()->GetSubMesh( thePrism.myBottom ); - if (( botSM->IsEmpty() ) && - ( ! botSM->GetAlgo() || - ! _gen->Compute( *botSM->GetFather(), botSM->GetSubShape(), /*shapeOnly=*/true ))) - return error( COMPERR_BAD_INPUT_MESH, - TCom( "No mesher defined to compute the base face #") - << shapeID( thePrism.myBottom )); + if ( !computeBase( thePrism )) + return false; // Make all side FACEs of thePrism meshed with quads if ( !computeWalls( thePrism )) @@ -1160,7 +1150,8 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism) // else if ( !trsf.empty() ) // bottomToTopTrsf = trsf.back(); - // To compute coordinates of a node inside a block, it is necessary to know + // To compute coordinates of a node inside a block using "block approach", + // it is necessary to know // 1. normalized parameters of the node by which // 2. coordinates of node projections on all block sub-shapes are computed @@ -1173,7 +1164,7 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism) // Projections on the top and bottom faces are taken from nodes existing // on these faces; find correspondence between bottom and top nodes - myUseBlock = false; + myUseBlock = false; // is set to true if projection is done using "block approach" myBotToColumnMap.clear(); if ( !assocOrProjBottom2Top( bottomToTopTrsf, thePrism ) ) // it also fills myBotToColumnMap return false; @@ -1181,38 +1172,62 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism) // Create nodes inside the block - // use transformation (issue 0020680, IPAL0052499) - StdMeshers_Sweeper sweeper; - double tol; - bool allowHighBndError; - if ( !myUseBlock ) { + // use transformation (issue 0020680, IPAL0052499) or a "straight line" approach + StdMeshers_Sweeper sweeper; + sweeper.myHelper = myHelper; + sweeper.myBotFace = thePrism.myBottom; + sweeper.myTopFace = thePrism.myTop; + // load boundary nodes into sweeper bool dummy; + std::set< const SMDS_MeshNode* > usedEndNodes; list< TopoDS_Edge >::const_iterator edge = thePrism.myBottomEdges.begin(); for ( ; edge != thePrism.myBottomEdges.end(); ++edge ) { int edgeID = meshDS->ShapeToIndex( *edge ); TParam2ColumnMap* u2col = const_cast ( myBlock.GetParam2ColumnMap( edgeID, dummy )); - TParam2ColumnMap::iterator u2colIt = u2col->begin(); - for ( ; u2colIt != u2col->end(); ++u2colIt ) + + TParam2ColumnMap::iterator u2colIt = u2col->begin(), u2colEnd = u2col->end(); + const SMDS_MeshNode* n0 = u2colIt->second[0]; + const SMDS_MeshNode* n1 = u2col->rbegin()->second[0]; + if ( !usedEndNodes.insert ( n0 ).second ) ++u2colIt; + if ( !usedEndNodes.insert ( n1 ).second ) --u2colEnd; + + for ( ; u2colIt != u2colEnd; ++u2colIt ) sweeper.myBndColumns.push_back( & u2colIt->second ); } - // load node columns inside the bottom face + // load node columns inside the bottom FACE TNode2ColumnMap::iterator bot_column = myBotToColumnMap.begin(); + sweeper.myIntColumns.reserve( myBotToColumnMap.size() ); for ( ; bot_column != myBotToColumnMap.end(); ++bot_column ) sweeper.myIntColumns.push_back( & bot_column->second ); - tol = getSweepTolerance( thePrism ); - allowHighBndError = !isSimpleBottom( thePrism ); - } + myHelper->SetElementsOnShape( true ); - if ( !myUseBlock && sweeper.ComputeNodes( *myHelper, tol, allowHighBndError )) - { + // If all "vertical" EDGEs are straight, then all nodes of an internal node column + // are located on a line connecting the top node and the bottom node. + bool isStrightColunm = allVerticalEdgesStraight( thePrism ); + if ( !isStrightColunm ) + { + double tol = getSweepTolerance( thePrism ); + bool allowHighBndError = !isSimpleBottom( thePrism ); + myUseBlock = !sweeper.ComputeNodesByTrsf( tol, allowHighBndError ); + } + else if ( sweeper.CheckSameZ() ) + { + myUseBlock = !sweeper.ComputeNodesOnStraightSameZ(); + } + else + { + myUseBlock = !sweeper.ComputeNodesOnStraight(); + } + myHelper->SetElementsOnShape( false ); } - else // use block approach + + if ( myUseBlock ) // use block approach { // loop on nodes inside the bottom face Prism_3D::TNode prevBNode; @@ -1361,6 +1376,75 @@ bool StdMeshers_Prism_3D::compute(const Prism_3D::TPrismTopo& thePrism) return true; } +//======================================================================= +//function : computeBase +//purpose : Compute the base face of a prism +//======================================================================= + +bool StdMeshers_Prism_3D::computeBase(const Prism_3D::TPrismTopo& thePrism) +{ + SMESH_Mesh* mesh = myHelper->GetMesh(); + SMESH_subMesh* botSM = mesh->GetSubMesh( thePrism.myBottom ); + if (( botSM->IsEmpty() ) && + ( ! botSM->GetAlgo() || + ! _gen->Compute( *botSM->GetFather(), botSM->GetSubShape(), /*shapeOnly=*/true ))) + { + // find any applicable algorithm assigned to any FACE of the main shape + std::vector< TopoDS_Shape > faces; + if ( myPrevBottomSM && + myPrevBottomSM->GetAlgo()->IsApplicableToShape( thePrism.myBottom, /*all=*/false )) + faces.push_back( myPrevBottomSM->GetSubShape() ); + + TopExp_Explorer faceIt( mesh->GetShapeToMesh(), TopAbs_FACE ); + for ( ; faceIt.More(); faceIt.Next() ) + faces.push_back( faceIt.Current() ); + + faces.push_back( TopoDS_Shape() ); // to try quadrangle algorithm + + SMESH_Algo* algo = 0; + for ( size_t i = 0; i < faces.size() && botSM->IsEmpty(); ++i ) + { + if ( faces[i].IsNull() ) algo = TQuadrangleAlgo::instance( this, myHelper ); + else algo = mesh->GetSubMesh( faces[i] )->GetAlgo(); + if ( algo && algo->IsApplicableToShape( thePrism.myBottom, /*all=*/false )) + { + // try to compute the bottom FACE + if ( algo->NeedDiscreteBoundary() ) + { + // compute sub-shapes + SMESH_subMeshIteratorPtr smIt = botSM->getDependsOnIterator(false,false); + bool subOK = true; + while ( smIt->more() && subOK ) + { + SMESH_subMesh* sub = smIt->next(); + sub->ComputeStateEngine( SMESH_subMesh::COMPUTE ); + subOK = sub->IsMeshComputed(); + } + if ( !subOK ) + continue; + } + try { + OCC_CATCH_SIGNALS; + algo->InitComputeError(); + algo->Compute( *mesh, botSM->GetSubShape() ); + } + catch (...) { + } + } + } + } + + if ( botSM->IsEmpty() ) + return error( COMPERR_BAD_INPUT_MESH, + TCom( "No mesher defined to compute the base face #") + << shapeID( thePrism.myBottom )); + + if ( botSM->GetAlgo() ) + myPrevBottomSM = botSM; + + return true; +} + //======================================================================= //function : computeWalls //purpose : Compute 2D mesh on walls FACEs of a prism @@ -1391,6 +1475,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) for ( ; quad != thePrism.myWallQuads[iW].end(); ++quad ) { StdMeshers_FaceSidePtr lftSide = (*quad)->side[ QUAD_LEFT_SIDE ]; + lftSide->Reverse(); // to go up for ( int i = 0; i < lftSide->NbEdges(); ++i ) { ++wgt[ iW ]; @@ -1418,6 +1503,11 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) for ( size_t iW = 0; iW != nbWalls; ++iW ) wgt2quad.insert( make_pair( wgt[ iW ], iW )); + // artificial quads to do outer <-> inner wall projection + std::map< int, FaceQuadStruct > iW2oiQuads; + std::map< int, FaceQuadStruct >::iterator w2oiq; + makeQuadsForOutInProjection( thePrism, wgt2quad, iW2oiQuads ); + // Project 'vertical' EDGEs, from left to right multimap< int, int >::reverse_iterator w2q = wgt2quad.rbegin(); for ( ; w2q != wgt2quad.rend(); ++w2q ) @@ -1434,10 +1524,25 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) if ( swapLeftRight ) std::swap( lftSide, rgtSide ); + bool isArtificialQuad = (( w2oiq = iW2oiQuads.find( iW )) != iW2oiQuads.end() ); + if ( isArtificialQuad ) + { + // reset sides to perform the outer <-> inner projection + FaceQuadStruct& oiQuad = w2oiq->second; + rgtSide = oiQuad.side[ QUAD_RIGHT_SIDE ]; + lftSide = oiQuad.side[ QUAD_LEFT_SIDE ]; + iW2oiQuads.erase( w2oiq ); + } + // assure that all the source (left) EDGEs are meshed int nbSrcSegments = 0; for ( int i = 0; i < lftSide->NbEdges(); ++i ) { + if ( isArtificialQuad ) + { + nbSrcSegments = lftSide->NbPoints()-1; + continue; + } const TopoDS_Edge& srcE = lftSide->Edge(i); SMESH_subMesh* srcSM = mesh->GetSubMesh( srcE ); if ( !srcSM->IsMeshComputed() ) { @@ -1513,7 +1618,7 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) const UVPtStructVec& srcNodeStr = lftSide->GetUVPtStruct(); if ( srcNodeStr.size() == 0 ) return toSM( error( TCom("Invalid node positions on edge #") << - shapeID( lftSide->Edge(0) ))); + lftSide->EdgeID(0) )); vector< SMDS_MeshNode* > newNodes( srcNodeStr.size() ); for ( int is2ndV = 0; is2ndV < 2; ++is2ndV ) { @@ -1521,12 +1626,12 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) TopoDS_Vertex v = myHelper->IthVertex( is2ndV, E ); mesh->GetSubMesh( v )->ComputeStateEngine( SMESH_subMesh::COMPUTE ); const SMDS_MeshNode* n = SMESH_Algo::VertexNode( v, meshDS ); - newNodes[ is2ndV ? 0 : newNodes.size()-1 ] = (SMDS_MeshNode*) n; + newNodes[ is2ndV ? newNodes.size()-1 : 0 ] = (SMDS_MeshNode*) n; } // compute nodes on target EDGEs DBGOUT( "COMPUTE V edge (proj) " << shapeID( lftSide->Edge(0))); - rgtSide->Reverse(); // direct it same as the lftSide + //rgtSide->Reverse(); // direct it same as the lftSide myHelper->SetElementsOnShape( false ); // myHelper holds the prism shape TopoDS_Edge tgtEdge; for ( size_t iN = 1; iN < srcNodeStr.size()-1; ++iN ) // add nodes @@ -1695,9 +1800,8 @@ bool StdMeshers_Prism_3D::computeWalls(const Prism_3D::TPrismTopo& thePrism) } //======================================================================= -/*! - * \brief Returns a source EDGE of propagation to a given EDGE - */ +//function : findPropagationSource +//purpose : Returns a source EDGE of propagation to a given EDGE //======================================================================= TopoDS_Edge StdMeshers_Prism_3D::findPropagationSource( const TopoDS_Edge& E ) @@ -1710,9 +1814,93 @@ TopoDS_Edge StdMeshers_Prism_3D::findPropagationSource( const TopoDS_Edge& E ) return TopoDS_Edge(); } +//======================================================================= +//function : makeQuadsForOutInProjection +//purpose : Create artificial wall quads for vertical projection between +// the outer and inner walls +//======================================================================= + +void StdMeshers_Prism_3D::makeQuadsForOutInProjection( const Prism_3D::TPrismTopo& thePrism, + multimap< int, int >& wgt2quad, + map< int, FaceQuadStruct >& iQ2oiQuads) +{ + if ( thePrism.NbWires() <= 1 ) + return; + + std::set< int > doneWires; // processed wires + + SMESH_Mesh* mesh = myHelper->GetMesh(); + const bool isForward = true; + const bool skipMedium = myHelper->GetIsQuadratic(); + + // make a source side for all projections + + multimap< int, int >::reverse_iterator w2q = wgt2quad.rbegin(); + const int iQuad = w2q->second; + const int iWire = getWireIndex( thePrism.myWallQuads[ iQuad ].front() ); + doneWires.insert( iWire ); + + UVPtStructVec srcNodes; + + Prism_3D::TQuadList::const_iterator quad = thePrism.myWallQuads[ iQuad ].begin(); + for ( ; quad != thePrism.myWallQuads[ iQuad ].end(); ++quad ) + { + StdMeshers_FaceSidePtr lftSide = (*quad)->side[ QUAD_LEFT_SIDE ]; + + // assure that all the source (left) EDGEs are meshed + for ( int i = 0; i < lftSide->NbEdges(); ++i ) + { + const TopoDS_Edge& srcE = lftSide->Edge(i); + SMESH_subMesh* srcSM = mesh->GetSubMesh( srcE ); + if ( !srcSM->IsMeshComputed() ) { + srcSM->ComputeSubMeshStateEngine( SMESH_subMesh::COMPUTE ); + srcSM->ComputeStateEngine ( SMESH_subMesh::COMPUTE ); + } + if ( !srcSM->IsMeshComputed() ) + return; + } + const UVPtStructVec& subNodes = lftSide->GetUVPtStruct(); + UVPtStructVec::const_iterator subBeg = subNodes.begin(), subEnd = subNodes.end(); + if ( !srcNodes.empty() ) ++subBeg; + srcNodes.insert( srcNodes.end(), subBeg, subEnd ); + } + StdMeshers_FaceSidePtr srcSide = StdMeshers_FaceSide::New( srcNodes ); + + // make the quads + + list< TopoDS_Edge > sideEdges; + TopoDS_Face face; + for ( ++w2q; w2q != wgt2quad.rend(); ++w2q ) + { + const int iQuad = w2q->second; + const Prism_3D::TQuadList& quads = thePrism.myWallQuads[ iQuad ]; + const int iWire = getWireIndex( quads.front() ); + if ( !doneWires.insert( iWire ).second ) + continue; + + sideEdges.clear(); + for ( quad = quads.begin(); quad != quads.end(); ++quad ) + { + StdMeshers_FaceSidePtr lftSide = (*quad)->side[ QUAD_LEFT_SIDE ]; + for ( int i = 0; i < lftSide->NbEdges(); ++i ) + sideEdges.push_back( lftSide->Edge( i )); + face = lftSide->Face(); + } + StdMeshers_FaceSidePtr tgtSide = + StdMeshers_FaceSide::New( face, sideEdges, mesh, isForward, skipMedium, myHelper ); + + FaceQuadStruct& newQuad = iQ2oiQuads[ iQuad ]; + newQuad.side.resize( 4 ); + newQuad.side[ QUAD_LEFT_SIDE ] = srcSide; + newQuad.side[ QUAD_RIGHT_SIDE ] = tgtSide; + + wgt2quad.insert( *w2q ); // to process this quad after processing the newQuad + } +} + //======================================================================= //function : Evaluate -//purpose : +//purpose : //======================================================================= bool StdMeshers_Prism_3D::Evaluate(SMESH_Mesh& theMesh, @@ -2274,7 +2462,7 @@ bool StdMeshers_Prism_3D::projectBottomToTop( const gp_Trsf & bottom // Check the projected mesh - if ( thePrism.myNbEdgesInWires.size() > 1 && // there are holes + if ( thePrism.NbWires() > 1 && // there are holes topHelper.IsDistorted2D( topSM, /*checkUV=*/false )) { SMESH_MeshEditor editor( topHelper.GetMesh() ); @@ -2400,7 +2588,7 @@ double StdMeshers_Prism_3D::getSweepTolerance( const Prism_3D::TPrismTopo& thePr bool StdMeshers_Prism_3D::isSimpleBottom( const Prism_3D::TPrismTopo& thePrism ) { - if ( thePrism.myBottomEdges.size() != 4 ) + if ( thePrism.myNbEdgesInWires.front() != 4 ) return false; // analyse angles between edges @@ -2432,6 +2620,43 @@ bool StdMeshers_Prism_3D::isSimpleBottom( const Prism_3D::TPrismTopo& thePrism ) return true; } +//======================================================================= +//function : allVerticalEdgesStraight +//purpose : Defines if all "vertical" EDGEs are straight +//======================================================================= + +bool StdMeshers_Prism_3D::allVerticalEdgesStraight( const Prism_3D::TPrismTopo& thePrism ) +{ + for ( size_t i = 0; i < thePrism.myWallQuads.size(); ++i ) + { + const Prism_3D::TQuadList& quads = thePrism.myWallQuads[i]; + Prism_3D::TQuadList::const_iterator quadIt = quads.begin(); + TopoDS_Edge prevQuadEdge; + for ( ; quadIt != quads.end(); ++quadIt ) + { + StdMeshers_FaceSidePtr rightSide = (*quadIt)->side[ QUAD_RIGHT_SIDE ]; + + if ( !prevQuadEdge.IsNull() && + !SMESH_Algo::IsContinuous( rightSide->Edge( 0 ), prevQuadEdge )) + return false; + + for ( int iE = 0; iE < rightSide->NbEdges(); ++iE ) + { + const TopoDS_Edge & rightE = rightSide->Edge( iE ); + if ( !SMESH_Algo::IsStraight( rightE, /*degenResult=*/true )) + return false; + + if ( iE > 0 && + !SMESH_Algo::IsContinuous( rightSide->Edge( iE-1 ), rightE )) + return false; + + prevQuadEdge = rightE; + } + } + } + return true; +} + //======================================================================= //function : project2dMesh //purpose : Project mesh faces from a source FACE of one prism (theSrcFace) @@ -3364,10 +3589,10 @@ bool StdMeshers_PrismAsBlock::Init(SMESH_MesherHelper* helper, if ( !myHelper->LoadNodeColumns( faceColumns, (*quad)->face, quadBot, meshDS )) return error(COMPERR_BAD_INPUT_MESH, TCom("Can't find regular quadrangle mesh ") << "on a side face #" << MeshDS()->ShapeToIndex( (*quad)->face )); - - if ( !faceColumns.empty() && (int)faceColumns.begin()->second.size() != VerticalSize() ) - return error(COMPERR_BAD_INPUT_MESH, "Different 'vertical' discretization"); } + if ( !faceColumns.empty() && (int)faceColumns.begin()->second.size() != VerticalSize() ) + return error(COMPERR_BAD_INPUT_MESH, "Different 'vertical' discretization"); + // edge columns int id = MeshDS()->ShapeToIndex( *edgeIt ); bool isForward = true; // meaningless for intenal wires @@ -4667,6 +4892,7 @@ bool StdMeshers_Sweeper::projectIntPoints(const vector< gp_XYZ >& fromBndPoin const vector< gp_XYZ >& toBndPoints, const vector< gp_XYZ >& fromIntPoints, vector< gp_XYZ >& toIntPoints, + const double r, NSProjUtils::TrsfFinder3D& trsf, vector< gp_XYZ > * bndError) { @@ -4691,54 +4917,34 @@ bool StdMeshers_Sweeper::projectIntPoints(const vector< gp_XYZ >& fromBndPoin (*bndError)[ iP ] = toBndPoints[ iP ] - fromTrsf; } } - return true; -} - -//================================================================================ -/*! - * \brief Add boundary error to ineternal points - */ -//================================================================================ -void StdMeshers_Sweeper::applyBoundaryError(const vector< gp_XYZ >& bndPoints, - const vector< gp_XYZ >& bndError1, - const vector< gp_XYZ >& bndError2, - const double r, - vector< gp_XYZ >& intPoints, - vector< double >& int2BndDist) -{ - // fix each internal point - const double eps = 1e-100; - for ( size_t iP = 0; iP < intPoints.size(); ++iP ) + // apply boundary error + if ( bndError && toIntPoints.size() == myTopBotTriangles.size() ) { - gp_XYZ & intPnt = intPoints[ iP ]; - - // compute distance from intPnt to each boundary node - double int2BndDistSum = 0; - for ( size_t iBnd = 0; iBnd < bndPoints.size(); ++iBnd ) + for ( size_t iP = 0; iP < toIntPoints.size(); ++iP ) { - int2BndDist[ iBnd ] = 1 / (( intPnt - bndPoints[ iBnd ]).SquareModulus() + eps ); - int2BndDistSum += int2BndDist[ iBnd ]; - } - - // apply bndError - for ( size_t iBnd = 0; iBnd < bndPoints.size(); ++iBnd ) - { - intPnt += bndError1[ iBnd ] * ( 1 - r ) * int2BndDist[ iBnd ] / int2BndDistSum; - intPnt += bndError2[ iBnd ] * r * int2BndDist[ iBnd ] / int2BndDistSum; + const TopBotTriangles& tbTrias = myTopBotTriangles[ iP ]; + for ( int i = 0; i < 3; ++i ) // boundary errors at 3 triangle nodes + { + toIntPoints[ iP ] += + ( (*bndError)[ tbTrias.myBotTriaNodes[i] ] * tbTrias.myBotBC[i] * ( 1 - r ) + + (*bndError)[ tbTrias.myTopTriaNodes[i] ] * tbTrias.myTopBC[i] * ( r )); + } } } + + return true; } //================================================================================ /*! - * \brief Creates internal nodes of the prism + * \brief Create internal nodes of the prism by computing an affine transformation + * from layer to layer */ //================================================================================ -bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, - const double tol, - const bool allowHighBndError) +bool StdMeshers_Sweeper::ComputeNodesByTrsf( const double tol, + const bool allowHighBndError) { const size_t zSize = myBndColumns[0]->size(); const size_t zSrc = 0, zTgt = zSize-1; @@ -4754,6 +4960,11 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, intPntsOfLayer[ zTgt ][ iP ] = intPoint( iP, zTgt ); } + // for each internal column find boundary nodes whose error to use for correction + prepareTopBotDelaunay(); + if ( !findDelaunayTriangles()) + return false; + // compute coordinates of internal nodes by projecting (transfroming) src and tgt // nodes towards the central layer @@ -4780,10 +4991,12 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, } if (! projectIntPoints( fromSrcBndPnts, toSrcBndPnts, intPntsOfLayer[ zS-1 ], intPntsOfLayer[ zS ], + zS / ( zSize - 1.), trsfOfLayer [ zS-1 ], & bndError[ zS-1 ])) return false; if (! projectIntPoints( fromTgtBndPnts, toTgtBndPnts, intPntsOfLayer[ zT+1 ], intPntsOfLayer[ zT ], + zT / ( zSize - 1.), trsfOfLayer [ zT+1 ], & bndError[ zT+1 ])) return false; @@ -4822,10 +5035,12 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, } if (! projectIntPoints( fromSrcBndPnts, toSrcBndPnts, intPntsOfLayer[ zS-1 ], centerSrcIntPnts, + zS / ( zSize - 1.), trsfOfLayer [ zS-1 ], & bndError[ zS-1 ])) return false; if (! projectIntPoints( fromTgtBndPnts, toTgtBndPnts, intPntsOfLayer[ zT+1 ], centerTgtIntPnts, + zT / ( zSize - 1.), trsfOfLayer [ zT+1 ], & bndError[ zT+1 ])) return false; @@ -4844,24 +5059,7 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, (intPntsOfLayer[ zS-1 ][ iP ] - centerTgtIntPnts[ iP ]).SquareModulus() < tol*tol; } - // Evaluate an error of boundary points - - bool bndErrorIsSmall = true; - for ( size_t iP = 0; ( iP < myBndColumns.size() && bndErrorIsSmall ); ++iP ) - { - double sumError = 0; - for ( size_t z = 1; z < zS; ++z ) // loop on layers - sumError += ( bndError[ z-1 ][ iP ].Modulus() + - bndError[ zSize-z ][ iP ].Modulus() ); - - bndErrorIsSmall = ( sumError < tol ); - } - - if ( !bndErrorIsSmall && !allowHighBndError ) - return false; - // compute final points on the central layer - std::vector< double > int2BndDist( myBndColumns.size() ); // work array of applyBoundaryError() double r = zS / ( zSize - 1.); if ( zS == zT ) { @@ -4870,11 +5068,6 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, intPntsOfLayer[ zS ][ iP ] = ( 1 - r ) * centerSrcIntPnts[ iP ] + r * centerTgtIntPnts[ iP ]; } - if ( !bndErrorIsSmall ) - { - applyBoundaryError( toSrcBndPnts, bndError[ zS-1 ], bndError[ zS+1 ], r, - intPntsOfLayer[ zS ], int2BndDist ); - } } else { @@ -4885,17 +5078,9 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, intPntsOfLayer[ zT ][ iP ] = r * intPntsOfLayer[ zT ][ iP ] + ( 1 - r ) * centerTgtIntPnts[ iP ]; } - if ( !bndErrorIsSmall ) - { - applyBoundaryError( toSrcBndPnts, bndError[ zS-1 ], bndError[ zS+1 ], r, - intPntsOfLayer[ zS ], int2BndDist ); - applyBoundaryError( toTgtBndPnts, bndError[ zT+1 ], bndError[ zT-1 ], r, - intPntsOfLayer[ zT ], int2BndDist ); - } } - centerIntErrorIsSmall = true; // 3D_mesh_Extrusion_00/A3 - bndErrorIsSmall = true; + //centerIntErrorIsSmall = true; // 3D_mesh_Extrusion_00/A3 if ( !centerIntErrorIsSmall ) { // Compensate the central error; continue adding projection @@ -4927,9 +5112,11 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, } projectIntPoints( fromSrcBndPnts, toSrcBndPnts, fromSrcIntPnts, toSrcIntPnts, + zS / ( zSize - 1.), trsfOfLayer[ zS+1 ], & srcBndError ); projectIntPoints( fromTgtBndPnts, toTgtBndPnts, fromTgtIntPnts, toTgtIntPnts, + zT / ( zSize - 1.), trsfOfLayer[ zT-1 ], & tgtBndError ); // if ( zS == zTgt - 1 ) @@ -4960,15 +5147,6 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, zTIntPnts[ iP ] = r * zTIntPnts[ iP ] + ( 1 - r ) * toTgtIntPnts[ iP ]; } - // compensate bnd error - if ( !bndErrorIsSmall ) - { - applyBoundaryError( toSrcBndPnts, srcBndError, bndError[ zS+1 ], r, - intPntsOfLayer[ zS ], int2BndDist ); - applyBoundaryError( toTgtBndPnts, tgtBndError, bndError[ zT-1 ], r, - intPntsOfLayer[ zT ], int2BndDist ); - } - fromSrcBndPnts.swap( toSrcBndPnts ); fromSrcIntPnts.swap( toSrcIntPnts ); fromTgtBndPnts.swap( toTgtBndPnts ); @@ -4976,27 +5154,8 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, } } // if ( !centerIntErrorIsSmall ) - else if ( !bndErrorIsSmall ) - { - zS = zSrc + 1; - zT = zTgt - 1; - for ( ; zS < zT; ++zS, --zT ) // vertical loop on layers - { - for ( size_t iP = 0; iP < myBndColumns.size(); ++iP ) - { - toSrcBndPnts[ iP ] = bndPoint( iP, zS ); - toTgtBndPnts[ iP ] = bndPoint( iP, zT ); - } - // compensate bnd error - applyBoundaryError( toSrcBndPnts, bndError[ zS-1 ], bndError[ zS-1 ], 0.5, - intPntsOfLayer[ zS ], int2BndDist ); - applyBoundaryError( toTgtBndPnts, bndError[ zT+1 ], bndError[ zT+1 ], 0.5, - intPntsOfLayer[ zT ], int2BndDist ); - } - } - // cout << "centerIntErrorIsSmall = " << centerIntErrorIsSmall<< endl; - // cout << "bndErrorIsSmall = " << bndErrorIsSmall<< endl; + //cout << "centerIntErrorIsSmall = " << centerIntErrorIsSmall<< endl; // Create nodes for ( size_t iP = 0; iP < myIntColumns.size(); ++iP ) @@ -5005,10 +5164,301 @@ bool StdMeshers_Sweeper::ComputeNodes( SMESH_MesherHelper& helper, for ( size_t z = zSrc + 1; z < zTgt; ++z ) // vertical loop on layers { const gp_XYZ & xyz = intPntsOfLayer[ z ][ iP ]; - if ( !( nodeCol[ z ] = helper.AddNode( xyz.X(), xyz.Y(), xyz.Z() ))) + if ( !( nodeCol[ z ] = myHelper->AddNode( xyz.X(), xyz.Y(), xyz.Z() ))) return false; } } return true; } + +//================================================================================ +/*! + * \brief Check if all nodes of each layers have same logical Z + */ +//================================================================================ + +bool StdMeshers_Sweeper::CheckSameZ() +{ + myZColumns.resize( myBndColumns.size() ); + fillZColumn( myZColumns[0], *myBndColumns[0] ); + + bool sameZ = true; + const double tol = 0.1 * 1./ myBndColumns[0]->size(); + + // check columns based on VERTEXes + + vector< int > vertexIndex; + vertexIndex.push_back( 0 ); + for ( size_t iC = 1; iC < myBndColumns.size() && sameZ; ++iC ) + { + if ( myBndColumns[iC]->front()->GetPosition()->GetDim() > 0 ) + continue; // not on VERTEX + + vertexIndex.push_back( iC ); + fillZColumn( myZColumns[iC], *myBndColumns[iC] ); + + for ( size_t iZ = 0; iZ < myZColumns[0].size() && sameZ; ++iZ ) + sameZ = ( Abs( myZColumns[0][iZ] - myZColumns[iC][iZ]) < tol ); + } + + // check columns based on EDGEs, one per EDGE + + for ( size_t i = 1; i < vertexIndex.size() && sameZ; ++i ) + { + if ( vertexIndex[i] - vertexIndex[i-1] < 2 ) + continue; + + int iC = ( vertexIndex[i] + vertexIndex[i-1] ) / 2; + fillZColumn( myZColumns[iC], *myBndColumns[iC] ); + + for ( size_t iZ = 0; iZ < myZColumns[0].size() && sameZ; ++iZ ) + sameZ = ( Abs( myZColumns[0][iZ] - myZColumns[iC][iZ]) < tol ); + } + + if ( sameZ ) + { + myZColumns.resize(1); + } + else + { + for ( size_t iC = 1; iC < myBndColumns.size(); ++iC ) + fillZColumn( myZColumns[iC], *myBndColumns[iC] ); + } + + return sameZ; +} + +//================================================================================ +/*! + * \brief Create internal nodes of the prism all located on straight lines with + * the same distribution along the lines. + */ +//================================================================================ + +bool StdMeshers_Sweeper::ComputeNodesOnStraightSameZ() +{ + TZColumn& z = myZColumns[0]; + + for ( size_t i = 0; i < myIntColumns.size(); ++i ) + { + TNodeColumn& nodes = *myIntColumns[i]; + SMESH_NodeXYZ n0( nodes[0] ), n1( nodes.back() ); + + for ( size_t iZ = 0; iZ < z.size(); ++iZ ) + { + gp_XYZ p = n0 * ( 1 - z[iZ] ) + n1 * z[iZ]; + nodes[ iZ+1 ] = myHelper->AddNode( p.X(), p.Y(), p.Z() ); + } + } + + return true; +} + +//================================================================================ +/*! + * \brief Create internal nodes of the prism all located on straight lines with + * different distributions along the lines. + */ +//================================================================================ + +bool StdMeshers_Sweeper::ComputeNodesOnStraight() +{ + prepareTopBotDelaunay(); + + const SMDS_MeshNode *botNode, *topNode; + const BRepMesh_Triangle *topTria; + double botBC[3], topBC[3]; // barycentric coordinates + int botTriaNodes[3], topTriaNodes[3]; + bool checkUV = true; + + int nbInternalNodes = myIntColumns.size(); + myBotDelaunay->InitTraversal( nbInternalNodes ); + + while (( botNode = myBotDelaunay->NextNode( botBC, botTriaNodes ))) + { + TNodeColumn* column = myIntColumns[ myNodeID2ColID( botNode->GetID() )]; + + // find a Delaunay triangle containing the topNode + topNode = column->back(); + gp_XY topUV = myHelper->GetNodeUV( myTopFace, topNode, NULL, &checkUV ); + // get a starting triangle basing on that top and bot boundary nodes have same index + topTria = myTopDelaunay->GetTriangleNear( botTriaNodes[0] ); + topTria = myTopDelaunay->FindTriangle( topUV, topTria, topBC, topTriaNodes ); + if ( !topTria ) + return false; + + // create nodes along a line + SMESH_NodeXYZ botP( botNode ), topP( topNode); + for ( size_t iZ = 0; iZ < myZColumns[0].size(); ++iZ ) + { + // use barycentric coordinates as weight of Z of boundary columns + double botZ = 0, topZ = 0; + for ( int i = 0; i < 3; ++i ) + { + botZ += botBC[i] * myZColumns[ botTriaNodes[i] ][ iZ ]; + topZ += topBC[i] * myZColumns[ topTriaNodes[i] ][ iZ ]; + } + double rZ = double( iZ + 1 ) / ( myZColumns[0].size() + 1 ); + double z = botZ * ( 1 - rZ ) + topZ * rZ; + gp_XYZ p = botP * ( 1 - z ) + topP * z; + (*column)[ iZ+1 ] = myHelper->AddNode( p.X(), p.Y(), p.Z() ); + } + } + + return myBotDelaunay->NbVisitedNodes() == nbInternalNodes; +} + +//================================================================================ +/*! + * \brief Compute Z of nodes of a straight column + */ +//================================================================================ + +void StdMeshers_Sweeper::fillZColumn( TZColumn& zColumn, + TNodeColumn& nodes ) +{ + if ( zColumn.size() == nodes.size() - 2 ) + return; + + gp_Pnt p0 = SMESH_NodeXYZ( nodes[0] ); + gp_Vec line( p0, SMESH_NodeXYZ( nodes.back() )); + double len2 = line.SquareMagnitude(); + + zColumn.resize( nodes.size() - 2 ); + for ( size_t i = 0; i < zColumn.size(); ++i ) + { + gp_Vec vec( p0, SMESH_NodeXYZ( nodes[ i+1] )); + zColumn[i] = ( line * vec ) / len2; // param [0,1] on the line + } +} + +//================================================================================ +/*! + * \brief Initialize *Delaunay members + */ +//================================================================================ + +void StdMeshers_Sweeper::prepareTopBotDelaunay() +{ + UVPtStructVec botUV( myBndColumns.size() ); + UVPtStructVec topUV( myBndColumns.size() ); + for ( size_t i = 0; i < myBndColumns.size(); ++i ) + { + TNodeColumn& nodes = *myBndColumns[i]; + botUV[i].node = nodes[0]; + botUV[i].SetUV( myHelper->GetNodeUV( myBotFace, nodes[0] )); + topUV[i].node = nodes.back(); + topUV[i].SetUV( myHelper->GetNodeUV( myTopFace, nodes.back() )); + botUV[i].node->setIsMarked( true ); + } + TopoDS_Edge dummyE; + SMESH_Mesh* mesh = myHelper->GetMesh(); + TSideVector botWires( 1, StdMeshers_FaceSide::New( botUV, myBotFace, dummyE, mesh )); + TSideVector topWires( 1, StdMeshers_FaceSide::New( topUV, myTopFace, dummyE, mesh )); + + // Delaunay mesh on the FACEs. + bool checkUV = false; + myBotDelaunay.reset( new NSProjUtils::Delaunay( botWires, checkUV )); + myTopDelaunay.reset( new NSProjUtils::Delaunay( topWires, checkUV )); + + if ( myHelper->GetIsQuadratic() ) + { + // mark all medium nodes of faces on botFace to avoid their treating + SMESHDS_SubMesh* smDS = myHelper->GetMeshDS()->MeshElements( myBotFace ); + SMDS_ElemIteratorPtr eIt = smDS->GetElements(); + while ( eIt->more() ) + { + const SMDS_MeshElement* e = eIt->next(); + for ( int i = e->NbCornerNodes(), nb = e->NbNodes(); i < nb; ++i ) + e->GetNode( i )->setIsMarked( true ); + } + } + + // map to get a node column by a bottom node + myNodeID2ColID.Clear(/*doReleaseMemory=*/false); + myNodeID2ColID.ReSize( myIntColumns.size() ); + + // un-mark nodes to treat (internal bottom nodes) to be returned by myBotDelaunay + for ( size_t i = 0; i < myIntColumns.size(); ++i ) + { + const SMDS_MeshNode* botNode = myIntColumns[i]->front(); + botNode->setIsMarked( false ); + myNodeID2ColID.Bind( botNode->GetID(), i ); + } +} + +//================================================================================ +/*! + * \brief For each internal node column, find Delaunay triangles including it + * and Barycentric Coordinates within the triangles. Fill in myTopBotTriangles + */ +//================================================================================ + +bool StdMeshers_Sweeper::findDelaunayTriangles() +{ + const SMDS_MeshNode *botNode, *topNode; + const BRepMesh_Triangle *topTria; + TopBotTriangles tbTrias; + bool checkUV = true; + + int nbInternalNodes = myIntColumns.size(); + myTopBotTriangles.resize( nbInternalNodes ); + + myBotDelaunay->InitTraversal( nbInternalNodes ); + + while (( botNode = myBotDelaunay->NextNode( tbTrias.myBotBC, tbTrias.myBotTriaNodes ))) + { + int colID = myNodeID2ColID( botNode->GetID() ); + TNodeColumn* column = myIntColumns[ colID ]; + + // find a Delaunay triangle containing the topNode + topNode = column->back(); + gp_XY topUV = myHelper->GetNodeUV( myTopFace, topNode, NULL, &checkUV ); + // get a starting triangle basing on that top and bot boundary nodes have same index + topTria = myTopDelaunay->GetTriangleNear( tbTrias.myBotTriaNodes[0] ); + topTria = myTopDelaunay->FindTriangle( topUV, topTria, + tbTrias.myTopBC, tbTrias.myTopTriaNodes ); + if ( !topTria ) + tbTrias.SetTopByBottom(); + + myTopBotTriangles[ colID ] = tbTrias; + } + + if ( myBotDelaunay->NbVisitedNodes() < nbInternalNodes ) + return false; + + myBotDelaunay.reset(); + myTopDelaunay.reset(); + myNodeID2ColID.Clear(); + + return true; +} + +//================================================================================ +/*! + * \brief Initialize fields + */ +//================================================================================ + +StdMeshers_Sweeper::TopBotTriangles::TopBotTriangles() +{ + myBotBC[0] = myBotBC[1] = myBotBC[2] = myTopBC[0] = myTopBC[1] = myTopBC[2] = 0.; + myBotTriaNodes[0] = myBotTriaNodes[1] = myBotTriaNodes[2] = 0; + myTopTriaNodes[0] = myTopTriaNodes[1] = myTopTriaNodes[2] = 0; +} + +//================================================================================ +/*! + * \brief Set top data equal to bottom data + */ +//================================================================================ + +void StdMeshers_Sweeper::TopBotTriangles::SetTopByBottom() +{ + for ( int i = 0; i < 3; ++i ) + { + myTopBC[i] = myBotBC[i]; + myTopTriaNodes[i] = myBotTriaNodes[0]; + } +} diff --git a/src/StdMeshers/StdMeshers_Prism_3D.hxx b/src/StdMeshers/StdMeshers_Prism_3D.hxx index ddfbff96a..62f0b9a14 100644 --- a/src/StdMeshers/StdMeshers_Prism_3D.hxx +++ b/src/StdMeshers/StdMeshers_Prism_3D.hxx @@ -29,21 +29,20 @@ #include "SMESH_StdMeshers.hxx" -#include "SMDS_MeshNode.hxx" -#include "SMDS_TypeOfPosition.hxx" #include "SMESHDS_Mesh.hxx" -#include "SMESH_Algo.hxx" #include "SMESH_Block.hxx" #include "SMESH_Comment.hxx" #include "SMESH_Mesh.hxx" #include "SMESH_MesherHelper.hxx" #include "SMESH_TypeDefs.hxx" #include "SMESH_subMesh.hxx" +#include "StdMeshers_ProjectionUtils.hxx" #include #include #include #include +#include #include #include #include @@ -56,10 +55,6 @@ namespace Prism_3D struct TNode; struct TPrismTopo; } -namespace StdMeshers_ProjectionUtils -{ - class TrsfFinder3D; -} class SMESHDS_SubMesh; class TopoDS_Edge; @@ -106,12 +101,14 @@ namespace Prism_3D TopoDS_Face myBottom; TopoDS_Face myTop; std::list< TopoDS_Edge > myBottomEdges; - std::vector< TQuadList> myWallQuads; // wall sides can be vertically composite + std::vector< TQuadList> myWallQuads; // wall sides can be vertically composite std::vector< int > myRightQuadIndex; // index of right neighbour wall quad std::list< int > myNbEdgesInWires; bool myNotQuadOnTop; + size_t NbWires() const { return myNbEdgesInWires.size(); } + void Clear(); void SetUpsideDown(); }; @@ -418,12 +415,22 @@ private: */ struct StdMeshers_Sweeper { + // input data + SMESH_MesherHelper* myHelper; + TopoDS_Face myBotFace; + TopoDS_Face myTopFace; std::vector< TNodeColumn* > myBndColumns; // boundary nodes + // output data std::vector< TNodeColumn* > myIntColumns; // internal nodes - bool ComputeNodes( SMESH_MesherHelper& helper, - const double tol, - const bool allowHighBndError ); + bool ComputeNodesByTrsf( const double tol, + const bool allowHighBndError ); + + bool CheckSameZ(); + + bool ComputeNodesOnStraightSameZ(); + + bool ComputeNodesOnStraight(); private: @@ -433,19 +440,36 @@ private: gp_XYZ intPoint( int iP, int z ) const { return SMESH_TNodeXYZ( (*myIntColumns[ iP ])[ z ]); } - static bool projectIntPoints(const std::vector< gp_XYZ >& fromBndPoints, - const std::vector< gp_XYZ >& toBndPoints, - const std::vector< gp_XYZ >& fromIntPoints, - std::vector< gp_XYZ >& toIntPoints, - StdMeshers_ProjectionUtils::TrsfFinder3D& trsf, - std::vector< gp_XYZ > * bndError); - - static void applyBoundaryError(const std::vector< gp_XYZ >& bndPoints, - const std::vector< gp_XYZ >& bndError1, - const std::vector< gp_XYZ >& bndError2, - const double r, - std::vector< gp_XYZ >& toIntPoints, - std::vector< double >& int2BndDist); + bool projectIntPoints(const std::vector< gp_XYZ >& fromBndPoints, + const std::vector< gp_XYZ >& toBndPoints, + const std::vector< gp_XYZ >& fromIntPoints, + std::vector< gp_XYZ >& toIntPoints, + const double r, + StdMeshers_ProjectionUtils::TrsfFinder3D& trsf, + std::vector< gp_XYZ > * bndError); + + typedef std::vector< double > TZColumn; + static void fillZColumn( TZColumn& zColumn, + TNodeColumn& nodes ); + + void prepareTopBotDelaunay(); + bool findDelaunayTriangles(); + + std::vector< TZColumn > myZColumns; // Z distribution of boundary nodes + + StdMeshers_ProjectionUtils::DelaunayPtr myTopDelaunay; + StdMeshers_ProjectionUtils::DelaunayPtr myBotDelaunay; + TColStd_DataMapOfIntegerInteger myNodeID2ColID; + + // top and bottom Delaulay triangles including an internal column + struct TopBotTriangles + { + double myBotBC[3], myTopBC[3]; // barycentric coordinates of a node within a triangle + int myBotTriaNodes[3], myTopTriaNodes[3]; // indices of boundary columns + TopBotTriangles(); + void SetTopByBottom(); + }; + std::vector< TopBotTriangles> myTopBotTriangles; }; // =============================================== @@ -486,6 +510,10 @@ public: SMESH_MesherHelper* helper); static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } private: @@ -508,11 +536,22 @@ public: */ bool compute(const Prism_3D::TPrismTopo& thePrism); + /*! + * \brief Compute the base face of a prism + */ + bool computeBase(const Prism_3D::TPrismTopo& thePrism); + /*! * \brief Compute 2D mesh on walls FACEs of a prism */ bool computeWalls(const Prism_3D::TPrismTopo& thePrism); + /*! + * \brief Create artificial wall quads for vertical projection between the outer and inner walls + */ + void makeQuadsForOutInProjection( const Prism_3D::TPrismTopo& thePrism, + std::multimap< int, int >& wgt2quad, + std::map< int, FaceQuadStruct >& iW2oiQuads); /*! * \brief Returns a source EDGE of propagation to a given EDGE */ @@ -545,6 +584,11 @@ public: */ bool isSimpleBottom( const Prism_3D::TPrismTopo& thePrism ); + /*! + * \brief Defines if all "vertical" EDGEs are straight + */ + bool allVerticalEdgesStraight( const Prism_3D::TPrismTopo& thePrism ); + /*! * \brief Project mesh faces from a source FACE of one prism to * a source FACE of another prism @@ -578,6 +622,7 @@ private: StdMeshers_PrismAsBlock myBlock; SMESH_MesherHelper* myHelper; + SMESH_subMesh* myPrevBottomSM; std::vector myShapeXYZ; // point on each sub-shape of the block diff --git a/src/StdMeshers/StdMeshers_ProjectionUtils.cxx b/src/StdMeshers/StdMeshers_ProjectionUtils.cxx index a0164ef01..11f2b5e97 100644 --- a/src/StdMeshers/StdMeshers_ProjectionUtils.cxx +++ b/src/StdMeshers/StdMeshers_ProjectionUtils.cxx @@ -28,6 +28,7 @@ #include "StdMeshers_ProjectionUtils.hxx" #include "SMDS_EdgePosition.hxx" +#include "SMDS_FacePosition.hxx" #include "SMESHDS_Mesh.hxx" #include "SMESH_Algo.hxx" #include "SMESH_Block.hxx" @@ -46,6 +47,7 @@ #include "utilities.h" #include +#include #include #include #include @@ -476,6 +478,26 @@ namespace { return !allBndEdges.empty(); } + /*! + * \brief Convertor used in Delaunay constructor + */ + struct SideVector2UVPtStructVec + { + std::vector< const UVPtStructVec* > _uvVecs; + + SideVector2UVPtStructVec( const TSideVector& wires ) + { + _uvVecs.resize( wires.size() ); + for ( size_t i = 0; i < wires.size(); ++i ) + _uvVecs[ i ] = & wires[i]->GetUVPtStruct(); + } + + operator const std::vector< const UVPtStructVec* > & () const + { + return _uvVecs; + } + }; + } // namespace //======================================================================= @@ -1098,6 +1120,11 @@ bool StdMeshers_ProjectionUtils::FindSubShapeAssociation(const TopoDS_Shape& the // take care of proper association of propagated edges bool same1 = edge1.IsSame( edges1.front() ); bool same2 = edge2.IsSame( edges2.front() ); + if ( !same1 && !same2 ) + { + same1 = ( edges1.back().Orientation() == edge1.Orientation() ); + same2 = ( edges2.back().Orientation() == edge2.Orientation() ); + } if ( same1 != same2 ) { reverseEdges(edges2, nbE); @@ -2787,4 +2814,146 @@ namespace StdMeshers_ProjectionUtils } return true; } -} + + //================================================================================ + /*! + * \brief triangulate the srcFace in 2D + * \param [in] srcWires - boundary of the src FACE + */ + //================================================================================ + + Morph::Morph(const TSideVector& srcWires): + _delaunay( srcWires, /*checkUV=*/true ) + { + _srcSubMesh = srcWires[0]->GetMesh()->GetSubMesh( srcWires[0]->Face() ); + } + + //================================================================================ + /*! + * \brief Move non-marked target nodes + * \param [in,out] tgtHelper - helper + * \param [in] tgtWires - boundary nodes of the target FACE; must be in the + * same order as the nodes in srcWires given in the constructor + * \param [in] src2tgtNodes - map of src -> tgt nodes + * \param [in] moveAll - to move all nodes; if \c false, move only non-marked nodes + * \return bool - Ok or not + */ + //================================================================================ + + bool Morph::Perform(SMESH_MesherHelper& tgtHelper, + const TSideVector& tgtWires, + Handle(ShapeAnalysis_Surface) tgtSurface, + const TNodeNodeMap& src2tgtNodes, + const bool moveAll) + { + // get tgt boundary points corresponding to src boundary nodes + size_t nbP = 0; + for ( size_t iW = 0; iW < tgtWires.size(); ++iW ) + nbP += tgtWires[iW]->NbPoints() - 1; // 1st and last points coincide + if ( nbP != _delaunay.GetBndNodes().size() ) + return false; + + std::vector< gp_XY > tgtUV( nbP ); + for ( size_t iW = 0, iP = 0; iW < tgtWires.size(); ++iW ) + { + const UVPtStructVec& tgtPnt = tgtWires[iW]->GetUVPtStruct(); + for ( int i = 0, nb = tgtPnt.size() - 1; i < nb; ++i, ++iP ) + { + tgtUV[ iP ] = tgtPnt[i].UV(); + } + } + + SMESHDS_Mesh* tgtMesh = tgtHelper.GetMeshDS(); + const SMDS_MeshNode *srcNode, *tgtNode; + + // un-mark internal src nodes in order iterate them using _delaunay + int nbSrcNodes = 0; + SMDS_NodeIteratorPtr nIt = _srcSubMesh->GetSubMeshDS()->GetNodes(); + if ( !nIt || !nIt->more() ) return true; + if ( moveAll ) + { + nbSrcNodes = _srcSubMesh->GetSubMeshDS()->NbNodes(); + while ( nIt->more() ) + nIt->next()->setIsMarked( false ); + } + else + { + while ( nIt->more() ) + nbSrcNodes += int( !nIt->next()->isMarked() ); + } + + // Move tgt nodes + + double bc[3]; // barycentric coordinates + int nodeIDs[3]; // nodes of a delaunay triangle + const SMDS_FacePosition* pos; + + _delaunay.InitTraversal( nbSrcNodes ); + + while (( srcNode = _delaunay.NextNode( bc, nodeIDs ))) + { + // compute new coordinates for a corresponding tgt node + gp_XY uvNew( 0., 0. ), nodeUV; + for ( int i = 0; i < 3; ++i ) + uvNew += bc[i] * tgtUV[ nodeIDs[i]]; + gp_Pnt xyz = tgtSurface->Value( uvNew ); + + // find and move tgt node + TNodeNodeMap::const_iterator n2n = src2tgtNodes.find( srcNode ); + if ( n2n == src2tgtNodes.end() ) continue; + tgtNode = n2n->second; + tgtMesh->MoveNode( tgtNode, xyz.X(), xyz.Y(), xyz.Z() ); + + if (( pos = dynamic_cast< const SMDS_FacePosition* >( tgtNode->GetPosition() ))) + const_cast( pos )->SetParameters( uvNew.X(), uvNew.Y() ); + + --nbSrcNodes; + } + + return nbSrcNodes == 0; + + } // Morph::Perform + + //======================================================================= + //function : Delaunay + //purpose : construct from face sides + //======================================================================= + + Delaunay::Delaunay( const TSideVector& wires, bool checkUV ): + SMESH_Delaunay( SideVector2UVPtStructVec( wires ), + TopoDS::Face( wires[0]->FaceHelper()->GetSubShape() ), + wires[0]->FaceHelper()->GetSubShapeID() ) + { + _wire = wires[0]; // keep a wire to assure _helper to keep alive + _helper = _wire->FaceHelper(); + _checkUVPtr = checkUV ? & _checkUV : 0; + } + + //======================================================================= + //function : Delaunay + //purpose : construct from UVPtStructVec's + //======================================================================= + + Delaunay::Delaunay( const std::vector< const UVPtStructVec* > & boundaryNodes, + SMESH_MesherHelper& faceHelper, + bool checkUV): + SMESH_Delaunay( boundaryNodes, + TopoDS::Face( faceHelper.GetSubShape() ), + faceHelper.GetSubShapeID() ) + { + _helper = & faceHelper; + _checkUVPtr = checkUV ? & _checkUV : 0; + } + + //======================================================================= + //function : getNodeUV + //purpose : + //======================================================================= + + gp_XY Delaunay::getNodeUV( const TopoDS_Face& face, const SMDS_MeshNode* node ) const + { + return _helper->GetNodeUV( face, node, 0, _checkUVPtr ); + } + + +} // namespace StdMeshers_ProjectionUtils diff --git a/src/StdMeshers/StdMeshers_ProjectionUtils.hxx b/src/StdMeshers/StdMeshers_ProjectionUtils.hxx index ced8febbd..aa3828443 100644 --- a/src/StdMeshers/StdMeshers_ProjectionUtils.hxx +++ b/src/StdMeshers/StdMeshers_ProjectionUtils.hxx @@ -31,10 +31,13 @@ #include "SMESH_StdMeshers.hxx" #include "SMDS_MeshElement.hxx" +#include "SMESH_Delaunay.hxx" +#include "StdMeshers_FaceSide.hxx" +#include #include -#include #include +#include #include #include #include @@ -51,6 +54,7 @@ class SMESH_Mesh; class SMESH_subMesh; class TopoDS_Shape; +//----------------------------------------------------------------------------------------- /*! * \brief Struct used instead of a sole TopTools_DataMapOfShapeShape to avoid * problems with bidirectional bindings @@ -91,6 +95,7 @@ namespace StdMeshers_ProjectionUtils TIDCompare> TNodeNodeMap; + //----------------------------------------------------------------------------------------- /*! * \brief Finds transformation between two sets of 2D points using * a least square approximation @@ -111,6 +116,7 @@ namespace StdMeshers_ProjectionUtils bool IsIdentity() const { return ( _trsf.Form() == gp_Identity ); } }; + //----------------------------------------------------------------------------------------- /*! * \brief Finds transformation between two sets of 3D points using * a least square approximation @@ -136,6 +142,56 @@ namespace StdMeshers_ProjectionUtils bool Invert(); }; + //----------------------------------------------------------------------------------------- + /*! + * \brief Create a Delaunay triangulation of nodes on a face boundary + * and provide exploration of nodes shared by elements lying on + * the face. For a returned node, also return a Delaunay triangle + * the node lies in and its Barycentric Coordinates within the triangle. + * Only non-marked nodes are visited. + * + * The main methods are defined in ../SMESHUtils/SMESH_Delaunay.hxx + */ + class Delaunay : public SMESH_Delaunay + { + public: + + Delaunay( const TSideVector& wires, bool checkUV = false ); + + Delaunay( const std::vector< const UVPtStructVec* > & boundaryNodes, + SMESH_MesherHelper& faceHelper, + bool checkUV = false); + + protected: + virtual gp_XY getNodeUV( const TopoDS_Face& face, const SMDS_MeshNode* node ) const; + + private: + SMESH_MesherHelper* _helper; + StdMeshers_FaceSidePtr _wire; + bool *_checkUVPtr, _checkUV; + }; + typedef boost::shared_ptr< Delaunay > DelaunayPtr; + + //----------------------------------------------------------------------------------------- + /*! + * \brief Morph mesh on the target FACE to lie within FACE boundary w/o distortion + */ + class Morph + { + Delaunay _delaunay; + SMESH_subMesh* _srcSubMesh; + public: + + Morph(const TSideVector& srcWires); + + bool Perform(SMESH_MesherHelper& tgtHelper, + const TSideVector& tgtWires, + Handle(ShapeAnalysis_Surface) tgtSurface, + const TNodeNodeMap& src2tgtNodes, + const bool moveAll); + }; + + //----------------------------------------------------------------------------------------- /*! * \brief Looks for association of all sub-shapes of two shapes * \param theShape1 - shape 1 diff --git a/src/StdMeshers/StdMeshers_Projection_2D.cxx b/src/StdMeshers/StdMeshers_Projection_2D.cxx index eeb20f009..64ce054c5 100644 --- a/src/StdMeshers/StdMeshers_Projection_2D.cxx +++ b/src/StdMeshers/StdMeshers_Projection_2D.cxx @@ -1133,7 +1133,7 @@ namespace { { SMESH_subMesh* faceSM = helper.GetMesh()->GetSubMesh( helper.GetSubShape() ); - if ( helper.IsDistorted2D( faceSM, /*checkUV=*/true )) + //if ( helper.IsDistorted2D( faceSM, /*checkUV=*/true )) { SMESH_MeshEditor editor( helper.GetMesh() ); SMESHDS_SubMesh* smDS = faceSM->GetSubMeshDS(); @@ -1179,240 +1179,6 @@ namespace { return true; } - typedef list< pair< const SMDS_MeshNode*, const BRepMesh_Triangle* > > TNodeTriaList; - - //================================================================================ - /*! - * \brief Add in-FACE nodes surrounding a given node to a queue - */ - //================================================================================ - - void addCloseNodes( const SMDS_MeshNode* srcNode, - const BRepMesh_Triangle* bmTria, - const int srcFaceID, - TNodeTriaList & noTriQueue ) - { - // find in-FACE nodes - SMDS_ElemIteratorPtr elems = srcNode->GetInverseElementIterator(SMDSAbs_Face); - while ( elems->more() ) - { - const SMDS_MeshElement* elem = elems->next(); - if ( elem->getshapeId() == srcFaceID ) - { - for ( int i = 0, nb = elem->NbNodes(); i < nb; ++i ) - { - const SMDS_MeshNode* n = elem->GetNode( i ); - if ( !n->isMarked() ) - noTriQueue.push_back( make_pair( n, bmTria )); - } - } - } - } - - //================================================================================ - /*! - * \brief Find a delauney triangle containing a given 2D point and return - * barycentric coordinates within the found triangle - */ - //================================================================================ - - const BRepMesh_Triangle* findTriangle( const gp_XY& uv, - const BRepMesh_Triangle* bmTria, - Handle(BRepMesh_DataStructureOfDelaun)& triaDS, - double bc[3] ) - { - int nodeIDs[3]; - gp_XY nodeUVs[3]; - int linkIDs[3]; - Standard_Boolean ori[3]; - - while ( bmTria ) - { - // check bmTria - - triaDS->ElementNodes( *bmTria, nodeIDs ); - nodeUVs[0] = triaDS->GetNode( nodeIDs[0] ).Coord(); - nodeUVs[1] = triaDS->GetNode( nodeIDs[1] ).Coord(); - nodeUVs[2] = triaDS->GetNode( nodeIDs[2] ).Coord(); - - SMESH_MeshAlgos::GetBarycentricCoords( uv, - nodeUVs[0], nodeUVs[1], nodeUVs[2], - bc[0], bc[1] ); - if ( bc[0] >= 0 && bc[1] >= 0 && bc[0] + bc[1] <= 1 ) - { - bc[2] = 1 - bc[0] - bc[1]; - return bmTria; - } - - // look for a neighbor triangle, which is adjacent to a link intersected - // by a segment( triangle center -> uv ) - - gp_XY gc = ( nodeUVs[0] + nodeUVs[1] + nodeUVs[2] ) / 3.; - gp_XY seg = uv - gc; - - bmTria->Edges( linkIDs, ori ); - int triaID = triaDS->IndexOf( *bmTria ); - bmTria = 0; - - for ( int i = 0; i < 3; ++i ) - { - const BRepMesh_PairOfIndex & triIDs = triaDS->ElementsConnectedTo( linkIDs[i] ); - if ( triIDs.Extent() < 2 ) - continue; // no neighbor triangle - - // check if a link intersects gc2uv - const BRepMesh_Edge & link = triaDS->GetLink( linkIDs[i] ); - const BRepMesh_Vertex & n1 = triaDS->GetNode( link.FirstNode() ); - const BRepMesh_Vertex & n2 = triaDS->GetNode( link.LastNode() ); - gp_XY uv1 = n1.Coord(); - gp_XY lin = n2.Coord() - uv1; // link direction - - double crossSegLin = seg ^ lin; - if ( Abs( crossSegLin ) < std::numeric_limits::min() ) - continue; // parallel - - double uSeg = ( uv1 - gc ) ^ lin / crossSegLin; - if ( 0. <= uSeg && uSeg <= 1. ) - { - bmTria = & triaDS->GetElement( triIDs.Index( 1 + ( triIDs.Index(1) == triaID ))); - break; - } - } - } - return bmTria; - } - - //================================================================================ - /*! - * \brief Morph mesh on the target face to lie within FACE boundary w/o distortion - * - * algo: - * - make a CDT on the src FACE - * - find a triangle containing a src node and get its barycentric coordinates - * - move the node to a point with the same barycentric coordinates in a corresponding - * tgt triangle - */ - //================================================================================ - - bool morph( SMESH_MesherHelper& tgtHelper, - const TopoDS_Face& tgtFace, - const TopoDS_Face& srcFace, - const TSideVector& tgtWires, - const TSideVector& srcWires, - const TAssocTool::TNodeNodeMap& src2tgtNodes ) - { - if ( srcWires.size() != tgtWires.size() ) return false; - if ( srcWires.size() == 1 ) return false; // tmp - - // count boundary points - int iP = 1, nbP = 0; - for ( size_t iW = 0; iW < srcWires.size(); ++iW ) - nbP += srcWires[iW]->NbPoints() - 1; // 1st and last points coincide - - // fill boundary points - BRepMesh::Array1OfVertexOfDelaun srcVert( 1, 1 + nbP ), tgtVert( 1, 1 + nbP ); - vector< const SMDS_MeshNode* > bndSrcNodes( nbP + 1 ); bndSrcNodes[0] = 0; - BRepMesh_Vertex v( 0, 0, BRepMesh_Frontier ); - for ( size_t iW = 0; iW < srcWires.size(); ++iW ) - { - const UVPtStructVec& srcPnt = srcWires[iW]->GetUVPtStruct(); - const UVPtStructVec& tgtPnt = tgtWires[iW]->GetUVPtStruct(); - if ( srcPnt.size() != tgtPnt.size() ) return false; - - for ( int i = 0, nb = srcPnt.size() - 1; i < nb; ++i, ++iP ) - { - bndSrcNodes[ iP ] = srcPnt[i].node; - srcPnt[i].node->setIsMarked( true ); - - v.ChangeCoord() = srcPnt[i].UV(); - srcVert( iP ) = v; - v.ChangeCoord() = tgtPnt[i].UV(); - tgtVert( iP ) = v; - } - } - // triangulate the srcFace in 2D - BRepMesh_Delaun delauney( srcVert ); - Handle(BRepMesh_DataStructureOfDelaun) triaDS = delauney.Result(); - - Handle(ShapeAnalysis_Surface) tgtSurface = tgtHelper.GetSurface( tgtFace ); - SMESHDS_Mesh* srcMesh = srcWires[0]->GetMesh()->GetMeshDS(); - SMESHDS_Mesh* tgtMesh = tgtHelper.GetMeshDS(); - const SMDS_MeshNode *srcNode, *tgtNode; - const BRepMesh_Triangle *bmTria; - - // un-mark internal src nodes; later we will mark moved nodes - SMDS_NodeIteratorPtr nIt = srcMesh->MeshElements( srcFace )->GetNodes(); - if ( !nIt || !nIt->more() ) return true; - while ( nIt->more() ) - ( srcNode = nIt->next() )->setIsMarked( false ); - - // initialize a queue of nodes with starting triangles - const int srcFaceID = srcNode->getshapeId(); - TNodeTriaList noTriQueue; - size_t iBndSrcN = 1; - for ( ; iBndSrcN < bndSrcNodes.size() && noTriQueue.empty(); ++iBndSrcN ) - { - // get a triangle - const BRepMesh::ListOfInteger & linkIds = triaDS->LinksConnectedTo( iBndSrcN ); - const BRepMesh_PairOfIndex & triaIds = triaDS->ElementsConnectedTo( linkIds.First() ); - const BRepMesh_Triangle& tria = triaDS->GetElement( triaIds.Index(1) ); - - addCloseNodes( bndSrcNodes[ iBndSrcN ], &tria, srcFaceID, noTriQueue ); - } - - // Move tgt nodes - - double bc[3]; // barycentric coordinates - int nodeIDs[3]; - bool checkUV = true; - const SMDS_FacePosition* pos; - - while ( !noTriQueue.empty() ) - { - srcNode = noTriQueue.front().first; - bmTria = noTriQueue.front().second; - noTriQueue.pop_front(); - if ( srcNode->isMarked() ) - continue; - srcNode->setIsMarked( true ); - - // find a delauney triangle containing the src node - gp_XY uv = tgtHelper.GetNodeUV( srcFace, srcNode, NULL, &checkUV ); - bmTria = findTriangle( uv, bmTria, triaDS, bc ); - if ( !bmTria ) - continue; - - // compute new coordinates for a corresponding tgt node - gp_XY uvNew( 0., 0. ), nodeUV; - triaDS->ElementNodes( *bmTria, nodeIDs ); - for ( int i = 0; i < 3; ++i ) - uvNew += bc[i] * tgtVert( nodeIDs[i]).Coord(); - gp_Pnt xyz = tgtSurface->Value( uvNew ); - - // find and move tgt node - TAssocTool::TNodeNodeMap::const_iterator n2n = src2tgtNodes.find( srcNode ); - if ( n2n == src2tgtNodes.end() ) continue; - tgtNode = n2n->second; - tgtMesh->MoveNode( tgtNode, xyz.X(), xyz.Y(), xyz.Z() ); - - if (( pos = dynamic_cast< const SMDS_FacePosition* >( tgtNode->GetPosition() ))) - const_cast( pos )->SetParameters( uvNew.X(), uvNew.Y() ); - - addCloseNodes( srcNode, bmTria, srcFaceID, noTriQueue ); - - // assure that all src nodes are visited - for ( ; iBndSrcN < bndSrcNodes.size() && noTriQueue.empty(); ++iBndSrcN ) - { - const BRepMesh::ListOfInteger & linkIds = triaDS->LinksConnectedTo( iBndSrcN ); - const BRepMesh_PairOfIndex & triaIds = triaDS->ElementsConnectedTo( linkIds.First() ); - const BRepMesh_Triangle& tria = triaDS->GetElement( triaIds.Index(1) ); - addCloseNodes( bndSrcNodes[ iBndSrcN ], &tria, srcFaceID, noTriQueue ); - } - } - - return true; - } - //======================================================================= /* * Set initial association of VERTEXes for the case of projection @@ -1614,11 +1380,12 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& TopoDS_Edge srcE1 = srcEdges.front(), tgtE1 = tgtEdges.front(); TopoDS_Shape srcE1bis = shape2ShapeMap( tgtE1 ); reverse = ( ! srcE1.IsSame( srcE1bis )); - if ( reverse && - //_sourceHypo->HasVertexAssociation() && + if ( ( reverse || srcE1.Orientation() != srcE1bis.Orientation() ) && nbEdgesInWires.front() > 2 && helper.IsRealSeam( tgtEdges.front() )) { + if ( srcE1.Orientation() != srcE1bis.Orientation() ) + reverse = true; // projection to a face with seam EDGE; pb is that GetOrderedEdges() // always puts a seam EDGE first (if possible) and as a result // we can't use only theReverse flag to correctly associate source @@ -1933,7 +1700,9 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& // ---------------------------------------------------------------- if ( helper.IsDistorted2D( tgtSubMesh, /*checkUV=*/false, &helper )) { - morph( helper, tgtFace, srcFace, tgtWires, srcWires, _src2tgtNodes ); + TAssocTool::Morph morph( srcWires ); + morph.Perform( helper, tgtWires, helper.GetSurface( tgtFace ), + _src2tgtNodes, /*moveAll=*/true ); if ( !fixDistortedFaces( helper, tgtWires )) return error("Invalid mesh generated"); diff --git a/src/StdMeshers/StdMeshers_Projection_3D.hxx b/src/StdMeshers/StdMeshers_Projection_3D.hxx index b1a05f9ac..008ca7429 100644 --- a/src/StdMeshers/StdMeshers_Projection_3D.hxx +++ b/src/StdMeshers/StdMeshers_Projection_3D.hxx @@ -57,6 +57,10 @@ public: */ virtual void SetEventListener(SMESH_subMesh* whenSetToSubMesh); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); protected: diff --git a/src/StdMeshers/StdMeshers_QuadFromMedialAxis_1D2D.hxx b/src/StdMeshers/StdMeshers_QuadFromMedialAxis_1D2D.hxx index 58a3364c7..43a62deba 100644 --- a/src/StdMeshers/StdMeshers_QuadFromMedialAxis_1D2D.hxx +++ b/src/StdMeshers/StdMeshers_QuadFromMedialAxis_1D2D.hxx @@ -53,6 +53,10 @@ class STDMESHERS_EXPORT StdMeshers_QuadFromMedialAxis_1D2D: public StdMeshers_Qu virtual void SetEventListener(SMESH_subMesh* subMesh); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); class Algo1D; diff --git a/src/StdMeshers/StdMeshers_Quadrangle_2D.hxx b/src/StdMeshers/StdMeshers_Quadrangle_2D.hxx index a1938bde9..b35a1425f 100644 --- a/src/StdMeshers/StdMeshers_Quadrangle_2D.hxx +++ b/src/StdMeshers/StdMeshers_Quadrangle_2D.hxx @@ -158,6 +158,10 @@ class STDMESHERS_EXPORT StdMeshers_Quadrangle_2D: public SMESH_2D_Algo const bool considerMesh = false, SMESH_MesherHelper* aFaceHelper = 0); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); protected: diff --git a/src/StdMeshers/StdMeshers_RadialPrism_3D.hxx b/src/StdMeshers/StdMeshers_RadialPrism_3D.hxx index 1f20b17fc..144905e81 100644 --- a/src/StdMeshers/StdMeshers_RadialPrism_3D.hxx +++ b/src/StdMeshers/StdMeshers_RadialPrism_3D.hxx @@ -55,6 +55,10 @@ public: virtual bool Evaluate(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape, MapShapeNbElems& aResMap); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); protected: diff --git a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.hxx b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.hxx index 7270cba40..a029757dd 100644 --- a/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.hxx +++ b/src/StdMeshers/StdMeshers_RadialQuadrangle_1D2D.hxx @@ -58,9 +58,13 @@ public: */ virtual void SubmeshRestored(SMESH_subMesh* subMesh); + virtual bool IsApplicableToShape(const TopoDS_Shape & shape, bool toCheckAll) const + { + return IsApplicable( shape, toCheckAll ); + } static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll); -protected: + protected: int computeLayerPositions(StdMeshers_FaceSidePtr linSide, std::vector< double >& positions, diff --git a/src/StdMeshers/StdMeshers_ViscousLayers.cxx b/src/StdMeshers/StdMeshers_ViscousLayers.cxx index a6530b924..8212b0a8b 100644 --- a/src/StdMeshers/StdMeshers_ViscousLayers.cxx +++ b/src/StdMeshers/StdMeshers_ViscousLayers.cxx @@ -95,7 +95,7 @@ #include #ifdef _DEBUG_ -//#define __myDEBUG +#define __myDEBUG //#define __NOT_INVALIDATE_BAD_SMOOTH //#define __NODES_AT_POS #endif @@ -427,20 +427,21 @@ namespace VISCOUS_3D enum EFlags { TO_SMOOTH = 0x0000001, MOVED = 0x0000002, // set by _neibors[i]->SetNewLength() - SMOOTHED = 0x0000004, // set by this->Smooth() + SMOOTHED = 0x0000004, // set by _LayerEdge::Smooth() DIFFICULT = 0x0000008, // near concave VERTEX ON_CONCAVE_FACE = 0x0000010, BLOCKED = 0x0000020, // not to inflate any more INTERSECTED = 0x0000040, // close intersection with a face found NORMAL_UPDATED = 0x0000080, - MARKED = 0x0000100, // local usage - MULTI_NORMAL = 0x0000200, // a normal is invisible by some of surrounding faces - NEAR_BOUNDARY = 0x0000400, // is near FACE boundary forcing smooth - SMOOTHED_C1 = 0x0000800, // is on _eosC1 - DISTORTED = 0x0001000, // was bad before smoothing - RISKY_SWOL = 0x0002000, // SWOL is parallel to a source FACE - SHRUNK = 0x0004000, // target node reached a tgt position while shrink() - UNUSED_FLAG = 0x0100000 // to add use flags after + UPD_NORMAL_CONV = 0x0000100, // to update normal on boundary of concave FACE + MARKED = 0x0000200, // local usage + MULTI_NORMAL = 0x0000400, // a normal is invisible by some of surrounding faces + NEAR_BOUNDARY = 0x0000800, // is near FACE boundary forcing smooth + SMOOTHED_C1 = 0x0001000, // is on _eosC1 + DISTORTED = 0x0002000, // was bad before smoothing + RISKY_SWOL = 0x0004000, // SWOL is parallel to a source FACE + SHRUNK = 0x0008000, // target node reached a tgt position while shrink() + UNUSED_FLAG = 0x0100000 // to add user flags after }; bool Is ( int flag ) const { return _flags & flag; } void Set ( int flag ) { _flags |= flag; } @@ -497,7 +498,7 @@ namespace VISCOUS_3D const gp_XYZ& PrevPos() const { return _pos[ _pos.size() - 2 ]; } gp_XYZ PrevCheckPos( _EdgesOnShape* eos=0 ) const; gp_Ax1 LastSegment(double& segLen, _EdgesOnShape& eos) const; - gp_XY LastUV( const TopoDS_Face& F, _EdgesOnShape& eos ) const; + gp_XY LastUV( const TopoDS_Face& F, _EdgesOnShape& eos, int which=-1 ) const; bool IsOnEdge() const { return _2neibors; } gp_XYZ Copy( _LayerEdge& other, _EdgesOnShape& eos, SMESH_MesherHelper& helper ); void SetCosin( double cosin ); @@ -664,6 +665,8 @@ namespace VISCOUS_3D _SolidData* _data; // parent SOLID + _LayerEdge* operator[](size_t i) const { return (_LayerEdge*) _edges[i]; } + size_t size() const { return _edges.size(); } TopAbs_ShapeEnum ShapeType() const { return _shape.IsNull() ? TopAbs_SHAPE : _shape.ShapeType(); } TopAbs_ShapeEnum SWOLType() const @@ -678,7 +681,7 @@ namespace VISCOUS_3D //-------------------------------------------------------------------------------- /*! - * \brief Convex FACE whose radius of curvature is less than the thickness of + * \brief Convex FACE whose radius of curvature is less than the thickness of * layers. It is used to detect distortion of prisms based on a convex * FACE and to update normals to enable further increasing the thickness */ @@ -692,7 +695,14 @@ namespace VISCOUS_3D // map a sub-shape to _SolidData::_edgesOnShape map< TGeomID, _EdgesOnShape* > _subIdToEOS; + bool _isTooCurved; bool _normalsFixed; + bool _normalsFixedOnBorders; // used in putOnOffsetSurface() + + double GetMaxCurvature( _SolidData& data, + _EdgesOnShape& eof, + BRepLProp_SLProps& surfProp, + SMESH_MesherHelper& helper); bool GetCenterOfCurvature( _LayerEdge* ledge, BRepLProp_SLProps& surfProp, @@ -757,8 +767,6 @@ namespace VISCOUS_3D int _nbShapesToSmooth; - //map< TGeomID,Handle(Geom_Curve)> _edge2curve; - vector< _CollisionEdges > _collisionEdges; set< TGeomID > _concaveFaces; @@ -870,6 +878,7 @@ namespace VISCOUS_3D const gp_XY& uvToFix, const double refSign ); }; + struct PyDump; //-------------------------------------------------------------------------------- /*! * \brief Builder of viscous layers @@ -942,8 +951,11 @@ namespace VISCOUS_3D void makeOffsetSurface( _EdgesOnShape& eos, SMESH_MesherHelper& ); void putOnOffsetSurface( _EdgesOnShape& eos, int infStep, vector< _EdgesOnShape* >& eosC1, - int smooStep=0, bool moveAll=false ); + int smooStep=0, int moveAll=false ); void findCollisionEdges( _SolidData& data, SMESH_MesherHelper& helper ); + void findEdgesToUpdateNormalNearConvexFace( _ConvexFace & convFace, + _SolidData& data, + SMESH_MesherHelper& helper ); void limitMaxLenByCurvature( _SolidData& data, SMESH_MesherHelper& helper ); void limitMaxLenByCurvature( _LayerEdge* e1, _LayerEdge* e2, _EdgesOnShape& eos1, _EdgesOnShape& eos2, @@ -987,6 +999,7 @@ namespace VISCOUS_3D TopTools_MapOfShape _shrinkedFaces; int _tmpFaceID; + PyDump* _pyDump; }; //-------------------------------------------------------------------------------- /*! @@ -1034,6 +1047,7 @@ namespace VISCOUS_3D size_t _iSeg[2]; // index of segment where extreme tgt node is projected _EdgesOnShape& _eos; double _curveLen; // length of the EDGE + std::pair _eToSmooth[2]; // indices of _LayerEdge's in _eos static Handle(Geom_Curve) CurveForSmooth( const TopoDS_Edge& E, _EdgesOnShape& eos, @@ -1047,31 +1061,24 @@ namespace VISCOUS_3D bool Perform(_SolidData& data, Handle(ShapeAnalysis_Surface)& surface, const TopoDS_Face& F, - SMESH_MesherHelper& helper ) - { - if ( _leParams.empty() || ( !isAnalytic() && _offPoints.empty() )) - prepare( data ); + SMESH_MesherHelper& helper ); - if ( isAnalytic() ) - return smoothAnalyticEdge( data, surface, F, helper ); - else - return smoothComplexEdge ( data, surface, F, helper ); - } void prepare(_SolidData& data ); + void findEdgesToSmooth(); + + bool isToSmooth( int iE ); + bool smoothAnalyticEdge( _SolidData& data, Handle(ShapeAnalysis_Surface)& surface, const TopoDS_Face& F, SMESH_MesherHelper& helper); - bool smoothComplexEdge( _SolidData& data, Handle(ShapeAnalysis_Surface)& surface, const TopoDS_Face& F, SMESH_MesherHelper& helper); - gp_XYZ getNormalNormal( const gp_XYZ & normal, const gp_XYZ& edgeDir); - _LayerEdge* getLEdgeOnV( bool is2nd ) { return _eos._edges[ is2nd ? _eos._edges.size()-1 : 0 ]->_2neibors->_edges[ is2nd ]; @@ -1663,14 +1670,15 @@ namespace VISCOUS_3D // HOWTO use: run python commands written in a console to see // construction steps of viscous layers #ifdef __myDEBUG - ofstream* py; - int theNbPyFunc; - struct PyDump { + ostream* py; + int theNbPyFunc; + struct PyDump + { PyDump(SMESH_Mesh& m) { int tag = 3 + m.GetId(); const char* fname = "/tmp/viscous.py"; cout << "execfile('"< ostream & operator<<( const T &anything ) { return *this ; } + }; + void Pause() { py = &_mystream; } + void Resume() { py = _pyStream; } + MyStream _mystream; + ostream* _pyStream; }; #define dumpFunction(f) { _dumpFunction(f, __LINE__);} #define dumpMove(n) { _dumpMove(n, __LINE__);} @@ -1710,7 +1726,7 @@ namespace VISCOUS_3D #else - struct PyDump { PyDump(SMESH_Mesh&) {} void Finish() {} }; + struct PyDump { PyDump(SMESH_Mesh&) {} void Finish() {} void Pause() {} void Resume() {} }; #define dumpFunction(f) f #define dumpMove(n) #define dumpMoveComm(n,txt) @@ -1852,6 +1868,7 @@ SMESH_ComputeErrorPtr _ViscousBuilder::Compute(SMESH_Mesh& theMesh, return SMESH_ComputeErrorPtr(); // everything already computed PyDump debugDump( theMesh ); + _pyDump = &debugDump; // TODO: ignore already computed SOLIDs if ( !findSolidsWithLayers()) @@ -2777,8 +2794,6 @@ void _ViscousBuilder::limitStepSizeByCurvature( _SolidData& data ) { SMESH_MesherHelper helper( *_mesh ); - const int nbTestPnt = 5; // on a FACE sub-shape - BRepLProp_SLProps surfProp( 2, 1e-6 ); data._convexFaces.clear(); @@ -2790,58 +2805,27 @@ void _ViscousBuilder::limitStepSizeByCurvature( _SolidData& data ) continue; TopoDS_Face F = TopoDS::Face( eof._shape ); - SMESH_subMesh * sm = eof._subMesh; const TGeomID faceID = eof._shapeID; BRepAdaptor_Surface surface( F, false ); surfProp.SetSurface( surface ); - bool isTooCurved = false; - _ConvexFace cnvFace; - const double oriFactor = ( F.Orientation() == TopAbs_REVERSED ? +1. : -1. ); - SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator(/*includeSelf=*/true); - while ( smIt->more() ) - { - sm = smIt->next(); - const TGeomID subID = sm->GetId(); - // find _LayerEdge's of a sub-shape - _EdgesOnShape* eos; - if (( eos = data.GetShapeEdges( subID ))) - cnvFace._subIdToEOS.insert( make_pair( subID, eos )); - else - continue; - // check concavity and curvature and limit data._stepSize - const double minCurvature = - 1. / ( eos->_hyp.GetTotalThickness() * ( 1 + theThickToIntersection )); - size_t iStep = Max( 1, eos->_edges.size() / nbTestPnt ); - for ( size_t i = 0; i < eos->_edges.size(); i += iStep ) - { - gp_XY uv = helper.GetNodeUV( F, eos->_edges[ i ]->_nodes[0] ); - surfProp.SetParameters( uv.X(), uv.Y() ); - if ( !surfProp.IsCurvatureDefined() ) - continue; - if ( surfProp.MaxCurvature() * oriFactor > minCurvature ) - { - limitStepSize( data, 0.9 / surfProp.MaxCurvature() * oriFactor ); - isTooCurved = true; - } - if ( surfProp.MinCurvature() * oriFactor > minCurvature ) - { - limitStepSize( data, 0.9 / surfProp.MinCurvature() * oriFactor ); - isTooCurved = true; - } - } - } // loop on sub-shapes of the FACE + cnvFace._face = F; + cnvFace._normalsFixed = false; + cnvFace._isTooCurved = false; - if ( !isTooCurved ) continue; + double maxCurvature = cnvFace.GetMaxCurvature( data, eof, surfProp, helper ); + if ( maxCurvature > 0 ) + { + limitStepSize( data, 0.9 / maxCurvature ); + findEdgesToUpdateNormalNearConvexFace( cnvFace, data, helper ); + } + if ( !cnvFace._isTooCurved ) continue; _ConvexFace & convFace = data._convexFaces.insert( make_pair( faceID, cnvFace )).first->second; - convFace._face = F; - convFace._normalsFixed = false; - // skip a closed surface (data._convexFaces is useful anyway) bool isClosedF = false; helper.SetSubShape( F ); @@ -2854,6 +2838,7 @@ void _ViscousBuilder::limitStepSizeByCurvature( _SolidData& data ) if ( isClosedF ) { // limit _LayerEdge::_maxLen on the FACE + const double oriFactor = ( F.Orientation() == TopAbs_REVERSED ? +1. : -1. ); const double minCurvature = 1. / ( eof._hyp.GetTotalThickness() * ( 1 + theThickToIntersection )); map< TGeomID, _EdgesOnShape* >::iterator id2eos = cnvFace._subIdToEOS.find( faceID ); @@ -2865,14 +2850,13 @@ void _ViscousBuilder::limitStepSizeByCurvature( _SolidData& data ) _LayerEdge* ledge = eos._edges[ i ]; gp_XY uv = helper.GetNodeUV( F, ledge->_nodes[0] ); surfProp.SetParameters( uv.X(), uv.Y() ); - if ( !surfProp.IsCurvatureDefined() ) - continue; - - if ( surfProp.MaxCurvature() * oriFactor > minCurvature ) - ledge->_maxLen = Min( ledge->_maxLen, 1. / surfProp.MaxCurvature() * oriFactor ); - - if ( surfProp.MinCurvature() * oriFactor > minCurvature ) - ledge->_maxLen = Min( ledge->_maxLen, 1. / surfProp.MinCurvature() * oriFactor ); + if ( surfProp.IsCurvatureDefined() ) + { + double curvature = Max( surfProp.MaxCurvature() * oriFactor, + surfProp.MinCurvature() * oriFactor ); + if ( curvature > minCurvature ) + ledge->_maxLen = Min( ledge->_maxLen, 1. / curvature ); + } } } continue; @@ -3049,8 +3033,8 @@ bool _ViscousBuilder::findShapesToSmooth( _SolidData& data ) { eos._edgeSmoother = new _Smoother1D( curve, eos ); - for ( size_t i = 0; i < eos._edges.size(); ++i ) - eos._edges[i]->Set( _LayerEdge::TO_SMOOTH ); + // for ( size_t i = 0; i < eos._edges.size(); ++i ) + // eos._edges[i]->Set( _LayerEdge::TO_SMOOTH ); } } } @@ -3429,11 +3413,12 @@ bool _ViscousBuilder::setEdgeData(_LayerEdge& edge, } // find _normal + bool fromVonF = false; if ( useGeometry ) { - bool fromVonF = ( eos.ShapeType() == TopAbs_VERTEX && - eos.SWOLType() == TopAbs_FACE && - totalNbFaces > 1 ); + fromVonF = ( eos.ShapeType() == TopAbs_VERTEX && + eos.SWOLType() == TopAbs_FACE && + totalNbFaces > 1 ); if ( onShrinkShape && !fromVonF ) // one of faces the node is on has no layers { @@ -3535,14 +3520,19 @@ bool _ViscousBuilder::setEdgeData(_LayerEdge& edge, break; } case TopAbs_VERTEX: { - //if ( eos.SWOLType() != TopAbs_FACE ) // else _cosin is set by getFaceDir() + if ( fromVonF ) + { + getFaceDir( TopoDS::Face( eos._sWOL ), TopoDS::Vertex( eos._shape ), + node, helper, normOK, &edge._cosin ); + } + else if ( eos.SWOLType() != TopAbs_FACE ) // else _cosin is set by getFaceDir() { TopoDS_Vertex V = TopoDS::Vertex( eos._shape ); gp_Vec inFaceDir = getFaceDir( F, V, node, helper, normOK ); double angle = inFaceDir.Angle( edge._normal ); // [0,PI] edge._cosin = Cos( angle ); if ( totalNbFaces > 2 || helper.IsSeamShape( node->getshapeId() )) - for ( int iF = totalNbFaces-2; iF >=0; --iF ) + for ( int iF = 1; iF < totalNbFaces; ++iF ) { F = face2Norm[ iF ].first; inFaceDir = getFaceDir( F, V, node, helper, normOK=true ); @@ -4428,11 +4418,14 @@ bool _ViscousBuilder::inflate(_SolidData& data) data._epsilon = data._stepSize * 1e-7; debugMsg( "-- geomSize = " << data._geomSize << ", stepSize = " << data._stepSize ); + _pyDump->Pause(); findCollisionEdges( data, helper ); limitMaxLenByCurvature( data, helper ); + _pyDump->Resume(); + // limit length of _LayerEdge's around MULTI_NORMAL _LayerEdge's for ( size_t i = 0; i < data._edgesOnShape.size(); ++i ) if ( data._edgesOnShape[i].ShapeType() == TopAbs_VERTEX && @@ -4959,7 +4952,7 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, // ignore intersection of a _LayerEdge based on a _ConvexFace with a face // lying on this _ConvexFace if ( _ConvexFace* convFace = data.GetConvexFace( intFace->getshapeId() )) - if ( convFace->_subIdToEOS.count ( eos._shapeID )) + if ( convFace->_isTooCurved && convFace->_subIdToEOS.count ( eos._shapeID )) continue; // ignore intersection of a _LayerEdge based on a FACE with an element on this FACE @@ -4971,15 +4964,16 @@ bool _ViscousBuilder::smoothAndCheck(_SolidData& data, if ( dist > 0 ) { bool toIgnore = false; - if ( eos._edges[i]->Is( _LayerEdge::TO_SMOOTH )) + if ( eos._toSmooth ) { const TopoDS_Shape& S = getMeshDS()->IndexToShape( intFace->getshapeId() ); if ( !S.IsNull() && S.ShapeType() == TopAbs_FACE ) { - TopExp_Explorer edge( eos._shape, TopAbs_EDGE ); - for ( ; !toIgnore && edge.More(); edge.Next() ) - // is adjacent - has a common EDGE - toIgnore = ( helper.IsSubShape( edge.Current(), S )); + TopExp_Explorer sub( eos._shape, + eos.ShapeType() == TopAbs_FACE ? TopAbs_EDGE : TopAbs_VERTEX ); + for ( ; !toIgnore && sub.More(); sub.Next() ) + // is adjacent - has a common EDGE or VERTEX + toIgnore = ( helper.IsSubShape( sub.Current(), S )); if ( toIgnore ) // check angle between normals { @@ -5264,7 +5258,7 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, int infStep, vector< _EdgesOnShape* >& eosC1, int smooStep, - bool moveAll ) + int moveAll ) { _EdgesOnShape * eof = & eos; if ( eos.ShapeType() != TopAbs_FACE ) // eos is a boundary of C1 FACE, look for the FACE eos @@ -5295,8 +5289,13 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, edge->Unset( _LayerEdge::MARKED ); if ( edge->Is( _LayerEdge::BLOCKED ) || !edge->_curvature ) continue; - if ( !moveAll && !edge->Is( _LayerEdge::MOVED )) + if ( moveAll == _LayerEdge::UPD_NORMAL_CONV ) + { + if ( !edge->Is( _LayerEdge::UPD_NORMAL_CONV )) continue; + } + else if ( !moveAll && !edge->Is( _LayerEdge::MOVED )) + continue; int nbBlockedAround = 0; for ( size_t iN = 0; iN < edge->_neibors.size(); ++iN ) @@ -5306,7 +5305,7 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, gp_Pnt tgtP = SMESH_TNodeXYZ( edge->_nodes.back() ); gp_Pnt2d uv = eof->_offsetSurf->NextValueOfUV( edge->_curvature->_uv, tgtP, preci ); - if ( eof->_offsetSurf->Gap() > edge->_len ) continue; // NextValueOfUV() bug + if ( eof->_offsetSurf->Gap() > edge->_len ) continue; // NextValueOfUV() bug edge->_curvature->_uv = uv; if ( eof->_offsetSurf->Gap() < 10 * preci ) continue; // same pos @@ -5325,9 +5324,15 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, edge->_pos.back() = newP; edge->Set( _LayerEdge::MARKED ); + if ( moveAll == _LayerEdge::UPD_NORMAL_CONV ) + { + edge->_normal = ( newP - prevP ).Normalized(); + } } } + + #ifdef _DEBUG_ // dumpMove() for debug size_t i = 0; @@ -5336,7 +5341,7 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, break; if ( i < eos._edges.size() ) { - dumpFunction(SMESH_Comment("putOnOffsetSurface_F") << eos._shapeID + dumpFunction(SMESH_Comment("putOnOffsetSurface_S") << eos._shapeID << "_InfStep" << infStep << "_" << smooStep ); for ( ; i < eos._edges.size(); ++i ) { @@ -5346,6 +5351,26 @@ void _ViscousBuilder::putOnOffsetSurface( _EdgesOnShape& eos, dumpFunctionEnd(); } #endif + + _ConvexFace* cnvFace; + if ( moveAll != _LayerEdge::UPD_NORMAL_CONV && + eos.ShapeType() == TopAbs_FACE && + (cnvFace = eos.GetData().GetConvexFace( eos._shapeID )) && + !cnvFace->_normalsFixedOnBorders ) + { + // put on the surface nodes built on FACE boundaries + SMESH_subMeshIteratorPtr smIt = eos._subMesh->getDependsOnIterator(/*includeSelf=*/false); + while ( smIt->more() ) + { + SMESH_subMesh* sm = smIt->next(); + _EdgesOnShape* subEOS = eos.GetData().GetShapeEdges( sm->GetId() ); + if ( !subEOS->_sWOL.IsNull() ) continue; + if ( std::find( eosC1.begin(), eosC1.end(), subEOS ) != eosC1.end() ) continue; + + putOnOffsetSurface( *subEOS, infStep, eosC1, smooStep, _LayerEdge::UPD_NORMAL_CONV ); + } + cnvFace->_normalsFixedOnBorders = true; + } } //================================================================================ @@ -5443,6 +5468,106 @@ Handle(Geom_Curve) _Smoother1D::CurveForSmooth( const TopoDS_Edge& E, return Handle(Geom_Curve)(); } +//================================================================================ +/*! + * \brief Smooth edges on EDGE + */ +//================================================================================ + +bool _Smoother1D::Perform(_SolidData& data, + Handle(ShapeAnalysis_Surface)& surface, + const TopoDS_Face& F, + SMESH_MesherHelper& helper ) +{ + if ( _leParams.empty() || ( !isAnalytic() && _offPoints.empty() )) + prepare( data ); + + findEdgesToSmooth(); + if ( isAnalytic() ) + return smoothAnalyticEdge( data, surface, F, helper ); + else + return smoothComplexEdge ( data, surface, F, helper ); +} + +//================================================================================ +/*! + * \brief Find edges to smooth + */ +//================================================================================ + +void _Smoother1D::findEdgesToSmooth() +{ + _LayerEdge* leOnV[2] = { getLEdgeOnV(0), getLEdgeOnV(1) }; + for ( int iEnd = 0; iEnd < 2; ++iEnd ) + if ( leOnV[iEnd]->Is( _LayerEdge::NORMAL_UPDATED )) + _leOnV[iEnd]._cosin = Abs( _edgeDir[iEnd].Normalized() * leOnV[iEnd]->_normal ); + + _eToSmooth[0].first = _eToSmooth[0].second = 0; + + for ( size_t i = 0; i < _eos.size(); ++i ) + { + if ( !_eos[i]->Is( _LayerEdge::TO_SMOOTH )) + { + if ( needSmoothing( _leOnV[0]._cosin, _eos[i]->_len, _curveLen * _leParams[i] ) || + isToSmooth( i )) + _eos[i]->Set( _LayerEdge::TO_SMOOTH ); + else + break; + } + _eToSmooth[0].second = i+1; + } + + _eToSmooth[1].first = _eToSmooth[1].second = _eos.size(); + + for ( int i = _eos.size() - 1; i >= _eToSmooth[0].second; --i ) + { + if ( !_eos[i]->Is( _LayerEdge::TO_SMOOTH )) + { + if ( needSmoothing( _leOnV[1]._cosin, _eos[i]->_len, _curveLen * ( 1.-_leParams[i] )) || + isToSmooth( i )) + _eos[i]->Set( _LayerEdge::TO_SMOOTH ); + else + break; + } + _eToSmooth[1].first = i; + } +} + +//================================================================================ +/*! + * \brief Check if iE-th _LayerEdge needs smoothing + */ +//================================================================================ + +bool _Smoother1D::isToSmooth( int iE ) +{ + SMESH_NodeXYZ pi( _eos[iE]->_nodes[0] ); + SMESH_NodeXYZ p0( _eos[iE]->_2neibors->srcNode(0) ); + SMESH_NodeXYZ p1( _eos[iE]->_2neibors->srcNode(1) ); + gp_XYZ seg0 = pi - p0; + gp_XYZ seg1 = p1 - pi; + gp_XYZ tangent = seg0 + seg1; + double tangentLen = tangent.Modulus(); + double segMinLen = Min( seg0.Modulus(), seg1.Modulus() ); + if ( tangentLen < std::numeric_limits::min() ) + return false; + tangent /= tangentLen; + + for ( size_t i = 0; i < _eos[iE]->_neibors.size(); ++i ) + { + _LayerEdge* ne = _eos[iE]->_neibors[i]; + if ( !ne->Is( _LayerEdge::TO_SMOOTH ) || + ne->_nodes.size() < 2 || + ne->_nodes[0]->GetPosition()->GetDim() != 2 ) + continue; + gp_XYZ edgeVec = SMESH_NodeXYZ( ne->_nodes.back() ) - SMESH_NodeXYZ( ne->_nodes[0] ); + double proj = edgeVec * tangent; + if ( needSmoothing( 1., proj, segMinLen )) + return true; + } + return false; +} + //================================================================================ /*! * \brief smooth _LayerEdge's on a staight EDGE or circular EDGE @@ -5456,80 +5581,100 @@ bool _Smoother1D::smoothAnalyticEdge( _SolidData& data, { if ( !isAnalytic() ) return false; - const size_t iFrom = 0, iTo = _eos._edges.size(); + size_t iFrom = 0, iTo = _eos._edges.size(); if ( _anaCurve->IsKind( STANDARD_TYPE( Geom_Line ))) { if ( F.IsNull() ) // 3D { - SMESH_TNodeXYZ p0 ( _eos._edges[iFrom]->_2neibors->tgtNode(0) ); - SMESH_TNodeXYZ p1 ( _eos._edges[iTo-1]->_2neibors->tgtNode(1) ); SMESH_TNodeXYZ pSrc0( _eos._edges[iFrom]->_2neibors->srcNode(0) ); SMESH_TNodeXYZ pSrc1( _eos._edges[iTo-1]->_2neibors->srcNode(1) ); - gp_XYZ newPos, lineDir = pSrc1 - pSrc0; - _LayerEdge* vLE0 = _eos._edges[iFrom]->_2neibors->_edges[0]; - _LayerEdge* vLE1 = _eos._edges[iTo-1]->_2neibors->_edges[1]; - bool shiftOnly = ( vLE0->Is( _LayerEdge::NORMAL_UPDATED ) || - vLE0->Is( _LayerEdge::BLOCKED ) || - vLE1->Is( _LayerEdge::NORMAL_UPDATED ) || - vLE1->Is( _LayerEdge::BLOCKED )); - for ( size_t i = iFrom; i < iTo; ++i ) - { - _LayerEdge* edge = _eos._edges[i]; - SMDS_MeshNode* tgtNode = const_cast( edge->_nodes.back() ); - newPos = p0 * ( 1. - _leParams[i] ) + p1 * _leParams[i]; - - if ( shiftOnly || edge->Is( _LayerEdge::NORMAL_UPDATED )) - { - gp_XYZ curPos = SMESH_TNodeXYZ ( tgtNode ); - double shift = ( lineDir * ( newPos - pSrc0 ) - - lineDir * ( curPos - pSrc0 )); - newPos = curPos + lineDir * shift / lineDir.SquareModulus(); - } - if ( edge->Is( _LayerEdge::BLOCKED )) + //const gp_XYZ lineDir = pSrc1 - pSrc0; + //_LayerEdge* vLE0 = getLEdgeOnV( 0 ); + //_LayerEdge* vLE1 = getLEdgeOnV( 1 ); + // bool shiftOnly = ( vLE0->Is( _LayerEdge::NORMAL_UPDATED ) || + // vLE0->Is( _LayerEdge::BLOCKED ) || + // vLE1->Is( _LayerEdge::NORMAL_UPDATED ) || + // vLE1->Is( _LayerEdge::BLOCKED )); + for ( int iEnd = 0; iEnd < 2; ++iEnd ) + { + iFrom = _eToSmooth[ iEnd ].first, iTo = _eToSmooth[ iEnd ].second; + if ( iFrom >= iTo ) continue; + SMESH_TNodeXYZ p0( _eos[iFrom]->_2neibors->tgtNode(0) ); + SMESH_TNodeXYZ p1( _eos[iTo-1]->_2neibors->tgtNode(1) ); + double param0 = ( iFrom == 0 ) ? 0. : _leParams[ iFrom-1 ]; + double param1 = _leParams[ iTo ]; + for ( size_t i = iFrom; i < iTo; ++i ) { - SMESH_TNodeXYZ pSrc( edge->_nodes[0] ); - double curThick = pSrc.SquareDistance( tgtNode ); - double newThink = ( pSrc - newPos ).SquareModulus(); - if ( newThink > curThick ) - continue; + _LayerEdge* edge = _eos[i]; + SMDS_MeshNode* tgtNode = const_cast( edge->_nodes.back() ); + double param = ( _leParams[i] - param0 ) / ( param1 - param0 ); + gp_XYZ newPos = p0 * ( 1. - param ) + p1 * param; + + // if ( shiftOnly || edge->Is( _LayerEdge::NORMAL_UPDATED )) + // { + // gp_XYZ curPos = SMESH_TNodeXYZ ( tgtNode ); + // double shift = ( lineDir * ( newPos - pSrc0 ) - + // lineDir * ( curPos - pSrc0 )); + // newPos = curPos + lineDir * shift / lineDir.SquareModulus(); + // } + if ( edge->Is( _LayerEdge::BLOCKED )) + { + SMESH_TNodeXYZ pSrc( edge->_nodes[0] ); + double curThick = pSrc.SquareDistance( tgtNode ); + double newThink = ( pSrc - newPos ).SquareModulus(); + if ( newThink > curThick ) + continue; + } + edge->_pos.back() = newPos; + tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() ); + dumpMove( tgtNode ); } - edge->_pos.back() = newPos; - tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() ); - dumpMove( tgtNode ); } } else // 2D { - _LayerEdge* e0 = getLEdgeOnV( 0 ); - _LayerEdge* e1 = getLEdgeOnV( 1 ); - gp_XY uv0 = e0->LastUV( F, *data.GetShapeEdges( e0 )); - gp_XY uv1 = e1->LastUV( F, *data.GetShapeEdges( e1 )); - if ( e0->_nodes.back() == e1->_nodes.back() ) // closed edge + _LayerEdge* eV0 = getLEdgeOnV( 0 ); + _LayerEdge* eV1 = getLEdgeOnV( 1 ); + gp_XY uvV0 = eV0->LastUV( F, *data.GetShapeEdges( eV0 )); + gp_XY uvV1 = eV1->LastUV( F, *data.GetShapeEdges( eV1 )); + if ( eV0->_nodes.back() == eV1->_nodes.back() ) // closed edge { int iPeriodic = helper.GetPeriodicIndex(); if ( iPeriodic == 1 || iPeriodic == 2 ) { - uv1.SetCoord( iPeriodic, helper.GetOtherParam( uv1.Coord( iPeriodic ))); - if ( uv0.Coord( iPeriodic ) > uv1.Coord( iPeriodic )) - std::swap( uv0, uv1 ); + uvV1.SetCoord( iPeriodic, helper.GetOtherParam( uvV1.Coord( iPeriodic ))); + if ( uvV0.Coord( iPeriodic ) > uvV1.Coord( iPeriodic )) + std::swap( uvV0, uvV1 ); } } - const gp_XY rangeUV = uv1 - uv0; - for ( size_t i = iFrom; i < iTo; ++i ) - { - if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue; - gp_XY newUV = uv0 + _leParams[i] * rangeUV; - _eos._edges[i]->_pos.back().SetCoord( newUV.X(), newUV.Y(), 0 ); - - gp_Pnt newPos = surface->Value( newUV.X(), newUV.Y() ); - SMDS_MeshNode* tgtNode = const_cast( _eos._edges[i]->_nodes.back() ); - tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() ); - dumpMove( tgtNode ); + for ( int iEnd = 0; iEnd < 2; ++iEnd ) + { + iFrom = _eToSmooth[ iEnd ].first, iTo = _eToSmooth[ iEnd ].second; + if ( iFrom >= iTo ) continue; + _LayerEdge* e0 = _eos[iFrom]->_2neibors->_edges[0]; + _LayerEdge* e1 = _eos[iTo-1]->_2neibors->_edges[1]; + gp_XY uv0 = ( e0 == eV0 ) ? uvV0 : e0->LastUV( F, _eos ); + gp_XY uv1 = ( e1 == eV1 ) ? uvV1 : e1->LastUV( F, _eos ); + double param0 = ( iFrom == 0 ) ? 0. : _leParams[ iFrom-1 ]; + double param1 = _leParams[ iTo ]; + const gp_XY rangeUV = uv1 - uv0; + for ( size_t i = iFrom; i < iTo; ++i ) + { + if ( _eos[i]->Is( _LayerEdge::BLOCKED )) continue; + double param = ( _leParams[i] - param0 ) / ( param1 - param0 ); + gp_XY newUV = uv0 + param * rangeUV; + _eos[i]->_pos.back().SetCoord( newUV.X(), newUV.Y(), 0 ); + + gp_Pnt newPos = surface->Value( newUV.X(), newUV.Y() ); + SMDS_MeshNode* tgtNode = const_cast( _eos[i]->_nodes.back() ); + tgtNode->setXYZ( newPos.X(), newPos.Y(), newPos.Z() ); + dumpMove( tgtNode ); - SMDS_FacePosition* pos = static_cast( tgtNode->GetPosition() ); - pos->SetUParameter( newUV.X() ); - pos->SetVParameter( newUV.Y() ); + SMDS_FacePosition* pos = static_cast( tgtNode->GetPosition() ); + pos->SetUParameter( newUV.X() ); + pos->SetVParameter( newUV.Y() ); + } } } return true; @@ -5568,9 +5713,10 @@ bool _Smoother1D::smoothAnalyticEdge( _SolidData& data, if ( uLast < 0 ) uLast += 2 * M_PI; - for ( size_t i = iFrom; i < iTo; ++i ) + for ( size_t i = 0; i < _eos.size(); ++i ) { - if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue; + if ( _eos[i]->Is( _LayerEdge::BLOCKED )) continue; + //if ( !_eos[i]->Is( _LayerEdge::TO_SMOOTH )) continue; double u = uLast * _leParams[i]; gp_Pnt p = ElCLib::Value( u, newCirc ); _eos._edges[i]->_pos.back() = p.XYZ(); @@ -5602,9 +5748,10 @@ bool _Smoother1D::smoothAnalyticEdge( _SolidData& data, gp_Ax2d axis( center, vec0 ); gp_Circ2d circ( axis, radius ); - for ( size_t i = iFrom; i < iTo; ++i ) + for ( size_t i = 0; i < _eos.size(); ++i ) { - if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue; + if ( _eos[i]->Is( _LayerEdge::BLOCKED )) continue; + //if ( !_eos[i]->Is( _LayerEdge::TO_SMOOTH )) continue; double newU = uLast * _leParams[i]; gp_Pnt2d newUV = ElCLib::Value( newU, circ ); _eos._edges[i]->_pos.back().SetCoord( newUV.X(), newUV.Y(), 0 ); @@ -5639,7 +5786,9 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, if ( _offPoints.empty() ) return false; + // ---------------------------------------------- // move _offPoints along normals of _LayerEdge's + // ---------------------------------------------- _LayerEdge* e[2] = { getLEdgeOnV(0), getLEdgeOnV(1) }; if ( e[0]->Is( _LayerEdge::NORMAL_UPDATED )) @@ -5683,7 +5832,9 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, } } + // ----------------------------------------------------------------- // project tgt nodes of extreme _LayerEdge's to the offset segments + // ----------------------------------------------------------------- if ( e[0]->Is( _LayerEdge::NORMAL_UPDATED )) _iSeg[0] = 0; if ( e[1]->Is( _LayerEdge::NORMAL_UPDATED )) _iSeg[1] = _offPoints.size()-2; @@ -5747,7 +5898,9 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, if ( e[1]->_normal * vDiv1.XYZ() < 0 ) e[1]->_len += d1; else e[1]->_len -= d1; + // --------------------------------------------------------------------------------- // compute normalized length of the offset segments located between the projections + // --------------------------------------------------------------------------------- size_t iSeg = 0, nbSeg = _iSeg[1] - _iSeg[0] + 1; vector< double > len( nbSeg + 1 ); @@ -5771,12 +5924,15 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, _offPoints[ _iSeg[0] ]._xyz = pExtreme[0].XYZ(); _offPoints[ _iSeg[1]+ 1]._xyz = pExtreme[1].XYZ(); + // ------------------------------------------------------------- // distribute tgt nodes of _LayerEdge's between the projections + // ------------------------------------------------------------- iSeg = 0; - for ( size_t i = 0; i < _eos._edges.size(); ++i ) + for ( size_t i = 0; i < _eos.size(); ++i ) { - if ( _eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue; + if ( _eos[i]->Is( _LayerEdge::BLOCKED )) continue; + //if ( !_eos[i]->Is( _LayerEdge::TO_SMOOTH )) continue; while ( iSeg+2 < len.size() && _leParams[i] > len[ iSeg+1 ] ) iSeg++; double r = ( _leParams[i] - len[ iSeg ]) / ( len[ iSeg+1 ] - len[ iSeg ]); @@ -5785,17 +5941,17 @@ bool _Smoother1D::smoothComplexEdge( _SolidData& data, if ( surface.IsNull() ) { - _eos._edges[i]->_pos.back() = p; + _eos[i]->_pos.back() = p; } else // project a new node position to a FACE { - gp_Pnt2d uv ( _eos._edges[i]->_pos.back().X(), _eos._edges[i]->_pos.back().Y() ); + gp_Pnt2d uv ( _eos[i]->_pos.back().X(), _eos[i]->_pos.back().Y() ); gp_Pnt2d uv2( surface->NextValueOfUV( uv, p, fTol )); p = surface->Value( uv2 ).XYZ(); - _eos._edges[i]->_pos.back().SetCoord( uv2.X(), uv2.Y(), 0 ); + _eos[i]->_pos.back().SetCoord( uv2.X(), uv2.Y(), 0 ); } - SMDS_MeshNode* tgtNode = const_cast( _eos._edges[i]->_nodes.back() ); + SMDS_MeshNode* tgtNode = const_cast( _eos[i]->_nodes.back() ); tgtNode->setXYZ( p.X(), p.Y(), p.Z() ); dumpMove( tgtNode ); } @@ -5836,8 +5992,20 @@ void _Smoother1D::prepare(_SolidData& data) double fullLen = _leParams.back() + pPrev.Distance( SMESH_TNodeXYZ( getLEdgeOnV(1)->_nodes[0])); for ( size_t i = 0; i < _leParams.size()-1; ++i ) _leParams[i] = _leParams[i+1] / fullLen; + _leParams.back() = 1.; } + _LayerEdge* leOnV[2] = { getLEdgeOnV(0), getLEdgeOnV(1) }; + + // get cosin to use in findEdgesToSmooth() + _edgeDir[0] = getEdgeDir( E, leOnV[0]->_nodes[0], data.GetHelper() ); + _edgeDir[1] = getEdgeDir( E, leOnV[1]->_nodes[0], data.GetHelper() ); + _leOnV[0]._cosin = Abs( leOnV[0]->_cosin ); + _leOnV[1]._cosin = Abs( leOnV[1]->_cosin ); + if ( _eos._sWOL.IsNull() ) // 3D + for ( int iEnd = 0; iEnd < 2; ++iEnd ) + _leOnV[iEnd]._cosin = Abs( _edgeDir[iEnd].Normalized() * leOnV[iEnd]->_normal ); + if ( isAnalytic() ) return; @@ -5864,8 +6032,6 @@ void _Smoother1D::prepare(_SolidData& data) _offPoints[i]._param = GCPnts_AbscissaPoint::Length( c3dAdaptor, u0, u ) / _curveLen; } - _LayerEdge* leOnV[2] = { getLEdgeOnV(0), getLEdgeOnV(1) }; - // set _2edges _offPoints [0]._2edges.set( &_leOnV[0], &_leOnV[0], 0.5, 0.5 ); _offPoints.back()._2edges.set( &_leOnV[1], &_leOnV[1], 0.5, 0.5 ); @@ -5903,9 +6069,6 @@ void _Smoother1D::prepare(_SolidData& data) int iLBO = _offPoints.size() - 2; // last but one - _edgeDir[0] = getEdgeDir( E, leOnV[0]->_nodes[0], data.GetHelper() ); - _edgeDir[1] = getEdgeDir( E, leOnV[1]->_nodes[0], data.GetHelper() ); - _leOnV[ 0 ]._normal = getNormalNormal( leOnV[0]->_normal, _edgeDir[0] ); _leOnV[ 1 ]._normal = getNormalNormal( leOnV[1]->_normal, _edgeDir[1] ); _leOnV[ 0 ]._len = 0; @@ -6490,6 +6653,80 @@ void _ViscousBuilder::findCollisionEdges( _SolidData& data, SMESH_MesherHelper& } } +//================================================================================ +/*! + * \brief Find _LayerEdge's located on boundary of a convex FACE whose normal + * will be updated at each inflation step + */ +//================================================================================ + +void _ViscousBuilder::findEdgesToUpdateNormalNearConvexFace( _ConvexFace & convFace, + _SolidData& data, + SMESH_MesherHelper& helper ) +{ + const TGeomID convFaceID = getMeshDS()->ShapeToIndex( convFace._face ); + const double preci = BRep_Tool::Tolerance( convFace._face ); + Handle(ShapeAnalysis_Surface) surface = helper.GetSurface( convFace._face ); + + bool edgesToUpdateFound = false; + + map< TGeomID, _EdgesOnShape* >::iterator id2eos = convFace._subIdToEOS.begin(); + for ( ; id2eos != convFace._subIdToEOS.end(); ++id2eos ) + { + _EdgesOnShape& eos = * id2eos->second; + if ( !eos._sWOL.IsNull() ) continue; + if ( !eos._hyp.ToSmooth() ) continue; + for ( size_t i = 0; i < eos._edges.size(); ++i ) + { + _LayerEdge* ledge = eos._edges[ i ]; + if ( ledge->Is( _LayerEdge::UPD_NORMAL_CONV )) continue; // already checked + if ( ledge->Is( _LayerEdge::MULTI_NORMAL )) continue; // not inflatable + + gp_XYZ tgtPos = ( SMESH_NodeXYZ( ledge->_nodes[0] ) + + ledge->_normal * ledge->_lenFactor * ledge->_maxLen ); + + // the normal must be updated if distance from tgtPos to surface is less than + // target thickness + + // find an initial UV for search of a projection of tgtPos to surface + const SMDS_MeshNode* nodeInFace = 0; + SMDS_ElemIteratorPtr fIt = ledge->_nodes[0]->GetInverseElementIterator(SMDSAbs_Face); + while ( fIt->more() && !nodeInFace ) + { + const SMDS_MeshElement* f = fIt->next(); + if ( convFaceID != f->getshapeId() ) continue; + + SMDS_ElemIteratorPtr nIt = f->nodesIterator(); + while ( nIt->more() && !nodeInFace ) + { + const SMDS_MeshElement* n = nIt->next(); + if ( n->getshapeId() == convFaceID ) + nodeInFace = static_cast< const SMDS_MeshNode* >( n ); + } + } + if ( !nodeInFace ) + continue; + gp_XY uv = helper.GetNodeUV( convFace._face, nodeInFace ); + + // projection + surface->NextValueOfUV( uv, tgtPos, preci ); + double dist = surface->Gap(); + if ( dist < 0.95 * ledge->_maxLen ) + { + ledge->Set( _LayerEdge::UPD_NORMAL_CONV ); + if ( !ledge->_curvature ) ledge->_curvature = new _Curvature; + ledge->_curvature->_uv.SetCoord( uv.X(), uv.Y() ); + edgesToUpdateFound = true; + } + } + } + + if ( !convFace._isTooCurved && edgesToUpdateFound ) + { + data._convexFaces.insert( make_pair( convFaceID, convFace )).first->second; + } +} + //================================================================================ /*! * \brief Modify normals of _LayerEdge's on EDGE's to avoid intersection with @@ -6942,6 +7179,8 @@ bool _ViscousBuilder::updateNormalsOfConvexFaces( _SolidData& data, for ( ; id2face != data._convexFaces.end(); ++id2face ) { _ConvexFace & convFace = (*id2face).second; + convFace._normalsFixedOnBorders = false; // to update at each inflation step + if ( convFace._normalsFixed ) continue; // already fixed if ( convFace.CheckPrisms() ) @@ -7304,6 +7543,59 @@ bool _ViscousBuilder::updateNormalsOfConvexFaces( _SolidData& data, return true; } +//================================================================================ +/*! + * \brief Return max curvature of a FACE + */ +//================================================================================ + +double _ConvexFace::GetMaxCurvature( _SolidData& data, + _EdgesOnShape& eof, + BRepLProp_SLProps& surfProp, + SMESH_MesherHelper& helper) +{ + double maxCurvature = 0; + + TopoDS_Face F = TopoDS::Face( eof._shape ); + + const int nbTestPnt = 5; + const double oriFactor = ( F.Orientation() == TopAbs_REVERSED ? +1. : -1. ); + SMESH_subMeshIteratorPtr smIt = eof._subMesh->getDependsOnIterator(/*includeSelf=*/true); + while ( smIt->more() ) + { + SMESH_subMesh* sm = smIt->next(); + const TGeomID subID = sm->GetId(); + + // find _LayerEdge's of a sub-shape + _EdgesOnShape* eos; + if (( eos = data.GetShapeEdges( subID ))) + this->_subIdToEOS.insert( make_pair( subID, eos )); + else + continue; + + // check concavity and curvature and limit data._stepSize + const double minCurvature = + 1. / ( eos->_hyp.GetTotalThickness() * ( 1 + theThickToIntersection )); + size_t iStep = Max( 1, eos->_edges.size() / nbTestPnt ); + for ( size_t i = 0; i < eos->_edges.size(); i += iStep ) + { + gp_XY uv = helper.GetNodeUV( F, eos->_edges[ i ]->_nodes[0] ); + surfProp.SetParameters( uv.X(), uv.Y() ); + if ( surfProp.IsCurvatureDefined() ) + { + double curvature = Max( surfProp.MaxCurvature() * oriFactor, + surfProp.MinCurvature() * oriFactor ); + maxCurvature = Max( maxCurvature, curvature ); + + if ( curvature > minCurvature ) + this->_isTooCurved = true; + } + } + } // loop on sub-shapes of the FACE + + return maxCurvature; +} + //================================================================================ /*! * \brief Finds a center of curvature of a surface at a _LayerEdge @@ -7586,13 +7878,14 @@ gp_Ax1 _LayerEdge::LastSegment(double& segLen, _EdgesOnShape& eos) const //================================================================================ /*! - * \brief Return the last position of the target node on a FACE. + * \brief Return the last (or \a which) position of the target node on a FACE. * \param [in] F - the FACE this _LayerEdge is inflated along + * \param [in] which - index of position * \return gp_XY - result UV */ //================================================================================ -gp_XY _LayerEdge::LastUV( const TopoDS_Face& F, _EdgesOnShape& eos ) const +gp_XY _LayerEdge::LastUV( const TopoDS_Face& F, _EdgesOnShape& eos, int which ) const { if ( F.IsSame( eos._sWOL )) // F is my FACE return gp_XY( _pos.back().X(), _pos.back().Y() ); @@ -7601,7 +7894,7 @@ gp_XY _LayerEdge::LastUV( const TopoDS_Face& F, _EdgesOnShape& eos ) const return gp_XY( 1e100, 1e100 ); // _sWOL is EDGE of F; _pos.back().X() is the last U on the EDGE - double f, l, u = _pos.back().X(); + double f, l, u = _pos[ which < 0 ? _pos.size()-1 : which ].X(); Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface( TopoDS::Edge(eos._sWOL), F, f,l); if ( !C2d.IsNull() && f <= u && u <= l ) return C2d->Value( u ).XY(); @@ -7679,7 +7972,7 @@ bool _LayerEdge::SegTriaInter( const gp_Ax1& lastSegment, * \param [in] eov - EOS of the VERTEX * \param [in] eos - EOS of the FACE * \param [in] step - inflation step - * \param [in,out] badSmooEdges - not untangled _LayerEdge's + * \param [in,out] badSmooEdges - tangled _LayerEdge's */ //================================================================================ @@ -9031,6 +9324,11 @@ void _LayerEdge::Block( _SolidData& data ) //if ( Is( BLOCKED )) return; Set( BLOCKED ); + SMESH_Comment msg( "#BLOCK shape="); + msg << data.GetShapeEdges( this )->_shapeID + << ", nodes " << _nodes[0]->GetID() << ", " << _nodes.back()->GetID(); + dumpCmd( msg + " -- BEGIN") + _maxLen = _len; std::queue<_LayerEdge*> queue; queue.push( this ); @@ -9056,6 +9354,8 @@ void _LayerEdge::Block( _SolidData& data ) //if ( edge->_nodes[0]->getshapeId() == neibor->_nodes[0]->getshapeId() ) viscous_layers_00/A3 { newMaxLen *= edge->_lenFactor / neibor->_lenFactor; + // newMaxLen *= Min( edge->_lenFactor / neibor->_lenFactor, + // neibor->_lenFactor / edge->_lenFactor ); } if ( neibor->_maxLen > newMaxLen ) { @@ -9073,6 +9373,7 @@ void _LayerEdge::Block( _SolidData& data ) } } } + dumpCmd( msg + " -- END") } //================================================================================ @@ -9227,6 +9528,7 @@ std::string _LayerEdge::DumpFlags() const case BLOCKED: dump << "BLOCKED"; break; case INTERSECTED: dump << "INTERSECTED"; break; case NORMAL_UPDATED: dump << "NORMAL_UPDATED"; break; + case UPD_NORMAL_CONV: dump << "UPD_NORMAL_CONV"; break; case MARKED: dump << "MARKED"; break; case MULTI_NORMAL: dump << "MULTI_NORMAL"; break; case NEAR_BOUNDARY: dump << "NEAR_BOUNDARY"; break; @@ -9262,7 +9564,7 @@ bool _ViscousBuilder::refine(_SolidData& data) double f,l, u = 0; gp_XY uv; vector< gp_XYZ > pos3D; - bool isOnEdge; + bool isOnEdge, isTooConvexFace = false; TGeomID prevBaseId = -1; TNode2Edge* n2eMap = 0; TNode2Edge::iterator n2e; @@ -9306,6 +9608,9 @@ bool _ViscousBuilder::refine(_SolidData& data) for ( size_t j = 0; j < eos._eosC1[i]->_edges.size(); ++j ) eos._eosC1[i]->_edges[j]->Set( _LayerEdge::SMOOTHED_C1 ); } + isTooConvexFace = false; + if ( _ConvexFace* cf = data.GetConvexFace( eos._shapeID )) + isTooConvexFace = cf->_isTooCurved; } vector< double > segLen; @@ -9321,8 +9626,8 @@ bool _ViscousBuilder::refine(_SolidData& data) if ( eos._sWOL.IsNull() ) { bool useNormal = true; - bool usePos = false; - bool smoothed = false; + bool usePos = false; + bool smoothed = false; double preci = 0.1 * edge._len; if ( eos._toSmooth && edge._pos.size() > 2 ) { @@ -9330,8 +9635,7 @@ bool _ViscousBuilder::refine(_SolidData& data) } if ( smoothed ) { - if ( !surface.IsNull() && - !data._convexFaces.count( eos._shapeID )) // edge smoothed on FACE + if ( !surface.IsNull() && !isTooConvexFace ) // edge smoothed on FACE { useNormal = usePos = false; gp_Pnt2d uv = helper.GetNodeUV( geomFace, edge._nodes[0] );