From b9d215ee2d335070116724cd20b2e89eba0a86ef Mon Sep 17 00:00:00 2001 From: Cesar CONOPOIMA Date: Thu, 21 Mar 2024 09:12:46 +0000 Subject: [PATCH] [bos #40649][CEA] Export Structured mesh in CGNS format for 2DQuadrangle and 3DHexahedron(i,j,k). --- idl/SMESH_Mesh.idl | 3 + src/DriverCGNS/CMakeLists.txt | 2 + src/DriverCGNS/DriverStructuredCGNS_Write.cxx | 439 +++++++++ src/DriverCGNS/DriverStructuredCGNS_Write.hxx | 70 ++ src/SMDS/SMDS_Mesh.hxx | 5 +- src/SMESH/SMESH_Mesh.cxx | 46 + src/SMESH/SMESH_Mesh.hxx | 3 + src/SMESHDS/CMakeLists.txt | 4 +- src/SMESHDS/SMESHDS_Mesh.cxx | 87 +- src/SMESHDS/SMESHDS_Mesh.hxx | 12 + src/SMESHGUI/SMESHGUI.cxx | 86 +- src/SMESHGUI/SMESH_msg_en.ts | 4 + src/SMESHGUI/SMESH_msg_fr.ts | 4 + src/SMESHGUI/SMESH_msg_ja.ts | 4 + src/SMESHUtils/CMakeLists.txt | 3 + src/SMESHUtils/SMESH_RegularGrid.cxx | 862 ++++++++++++++++++ src/SMESHUtils/SMESH_RegularGrid.hxx | 363 ++++++++ src/SMESHUtils/SMESH_RegularGridTemplate.hxx | 63 ++ src/SMESH_I/SMESH_Mesh_i.cxx | 46 +- src/SMESH_I/SMESH_Mesh_i.hxx | 7 +- src/SMESH_SWIG/smeshBuilder.py | 20 + src/StdMeshers/StdMeshers_Hexa_3D.cxx | 18 +- src/StdMeshers/StdMeshers_Quadrangle_2D.cxx | 10 + test/CMakeLists.txt | 19 + test/CTestTestfileInstall.cmake | 7 + test/SMESH_RegularGridTest.cxx | 326 +++++++ test/tests.set | 5 + 27 files changed, 2485 insertions(+), 33 deletions(-) create mode 100644 src/DriverCGNS/DriverStructuredCGNS_Write.cxx create mode 100644 src/DriverCGNS/DriverStructuredCGNS_Write.hxx create mode 100644 src/SMESHUtils/SMESH_RegularGrid.cxx create mode 100644 src/SMESHUtils/SMESH_RegularGrid.hxx create mode 100644 src/SMESHUtils/SMESH_RegularGridTemplate.hxx create mode 100644 test/SMESH_RegularGridTest.cxx diff --git a/idl/SMESH_Mesh.idl b/idl/SMESH_Mesh.idl index 9a888182c..6535888a1 100644 --- a/idl/SMESH_Mesh.idl +++ b/idl/SMESH_Mesh.idl @@ -756,6 +756,9 @@ module SMESH in string file, in boolean overwrite, in boolean groupElemsByType) raises (SALOME::SALOME_Exception); + void ExportStructuredCGNS( in SMESH_IDSource meshPart, + in string file, + in boolean overwrite) raises (SALOME::SALOME_Exception); void ExportGMF( in SMESH_IDSource meshPart, in string file, in boolean withRequiredGroups) raises (SALOME::SALOME_Exception); diff --git a/src/DriverCGNS/CMakeLists.txt b/src/DriverCGNS/CMakeLists.txt index 0380176e4..d1c716355 100644 --- a/src/DriverCGNS/CMakeLists.txt +++ b/src/DriverCGNS/CMakeLists.txt @@ -50,6 +50,7 @@ SET(_link_LIBRARIES SET(MeshDriverCGNS_HEADERS DriverCGNS_Read.hxx DriverCGNS_Write.hxx + DriverStructuredCGNS_Write.hxx SMESH_DriverCGNS.hxx ) @@ -59,6 +60,7 @@ SET(MeshDriverCGNS_HEADERS SET(MeshDriverCGNS_SOURCES DriverCGNS_Read.cxx DriverCGNS_Write.cxx + DriverStructuredCGNS_Write.cxx ) # --- rules --- diff --git a/src/DriverCGNS/DriverStructuredCGNS_Write.cxx b/src/DriverCGNS/DriverStructuredCGNS_Write.cxx new file mode 100644 index 000000000..80e0b4006 --- /dev/null +++ b/src/DriverCGNS/DriverStructuredCGNS_Write.cxx @@ -0,0 +1,439 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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 : DriverStructuredCGNS_Write.cxx +// Created : Tuesday March 19 2024 +// Author : Cesar Conopoima (cce) + +#include "DriverStructuredCGNS_Write.hxx" + +#include "SMDS_IteratorOnIterators.hxx" +#include "SMDS_MeshNode.hxx" +#include "SMDS_VolumeTool.hxx" +#include "SMESHDS_GroupBase.hxx" +#include "SMESHDS_GroupOnGeom.hxx" +#include "SMESHDS_Mesh.hxx" +#include "SMESH_Block.hxx" +#include "SMESH_Comment.hxx" + +#include +#include +#include +#include +#include +#include + +//CAS +#include + +// CGNS +#include + +#if CGNS_VERSION < 3100 +# define cgsize_t int +#endif + + +std::string DriverStructuredCGNS_Write::GetGroupName( const int shapeToIndex, int dim ) +{ + auto groupBase = myMesh->GetGroups(); + for( auto grp : groupBase ) + { + /*associate the group based on geometry*/ + if ( dynamic_cast( grp )) + { + SMESHDS_GroupOnGeom * geomGrp = dynamic_cast( grp ); + TopoDS_Shape shape = geomGrp->GetShape(); + if( myMesh->ShapeToIndex(shape) == shapeToIndex ) + return std::string(geomGrp->GetStoreName()); + + if ( dim==3 && (shape.ShapeType() == TopAbs_COMPSOLID || shape.ShapeType() == TopAbs_COMPOUND )) + { + for ( TopExp_Explorer fEx( shape, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + if( myMesh->ShapeToIndex(fEx.Current()) == shapeToIndex ) + return std::string(geomGrp->GetStoreName()); + } + } + + if ( dim==2 && (shape.ShapeType() == TopAbs_SHELL || shape.ShapeType() == TopAbs_COMPOUND) ) + { + for ( TopExp_Explorer fEx( shape, TopAbs_FACE, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + if( myMesh->ShapeToIndex(fEx.Current()) == shapeToIndex ) + return std::string(geomGrp->GetStoreName()); + } + } + } + } + return std::string(""); +} + +void DriverStructuredCGNS_Write::CheckForGroupNameOnFaceInterfaces( const SMESHUtils::SMESH_RegularGrid* grid, std::vector& boundaryNames ) +{ + auto groupBase = myMesh->GetGroups(); + for( auto grp : groupBase ) + { + if ( dynamic_cast( grp ) ) + { + SMESHDS_GroupOnGeom * geomGrp = dynamic_cast( grp ); + TopoDS_Shape shape = geomGrp->GetShape(); + if( shape.ShapeType() == TopAbs_FACE ) + { + SMESHUtils::SMESH_RegularGrid::FaceType f = grid->getFaceTypeByGeomFace( shape ); + if ( f != SMESHUtils::SMESH_RegularGrid::B_NONE ) + boundaryNames[ int(f) ] = geomGrp->GetStoreName(); + } + + if ( shape.ShapeType() == TopAbs_SHELL || shape.ShapeType() == TopAbs_COMPOUND ) + { + for ( TopExp_Explorer fEx( shape, TopAbs_FACE, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + SMESHUtils::SMESH_RegularGrid::FaceType f = grid->getFaceTypeByGeomFace( fEx.Current() ); + if ( f != SMESHUtils::SMESH_RegularGrid::B_NONE ) + boundaryNames[ int(f) ] = geomGrp->GetStoreName(); + } + } + } + } +} + +void DriverStructuredCGNS_Write::CheckForGroupNameOnEdgeInterfaces( const SMESHUtils::SMESH_RegularGrid* grid, std::vector& boundaryNames ) +{ + auto groupBase = myMesh->GetGroups(); + for( auto grp : groupBase ) + { + if ( dynamic_cast( grp ) ) + { + SMESHDS_GroupOnGeom * geomGrp = dynamic_cast( grp ); + TopoDS_Shape shape = geomGrp->GetShape(); + if( shape.ShapeType() == TopAbs_EDGE ) + { + SMESHUtils::SMESH_RegularGrid::EdgeType e = grid->getEdgeTypeByGeomEdge( shape ); + if ( e != SMESHUtils::SMESH_RegularGrid::NONE ) + boundaryNames[ int(e) ] = geomGrp->GetStoreName(); + } + + if ( shape.ShapeType() == TopAbs_COMPOUND ) + { + for ( TopExp_Explorer fEx( shape, TopAbs_EDGE ); fEx.More(); fEx.Next() ) + { + SMESHUtils::SMESH_RegularGrid::EdgeType e = grid->getEdgeTypeByGeomEdge( fEx.Current() ); + if ( e != SMESHUtils::SMESH_RegularGrid::NONE ) + boundaryNames[ int(e) ] = geomGrp->GetStoreName(); + } + } + } + } +} + +//================================================================================ +/*! + * \brief Write the mesh into the CGNS file + */ +//================================================================================ + +Driver_Mesh::Status DriverStructuredCGNS_Write::Perform() +{ + myErrorMessages.clear(); + + if ( !myMesh || myMesh->GetMeshInfo().NbElements() < 1 ) + return addMessage( !myMesh ? "NULL mesh" : "Empty mesh (no elements)", /*fatal = */true ); + + if ( Driver_Mesh::IsMeshTooLarge< cgsize_t >( myMesh, /*checkIDs =*/ false)) + return DRS_TOO_LARGE_MESH; + + // open the file + if ( cg_open(myFile.c_str(), CG_MODE_MODIFY, &_fn) != CG_OK && + cg_open(myFile.c_str(), CG_MODE_WRITE, &_fn) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // create a Base + // -------------- + + const int spaceDim = 3; + int meshDim = 1; + if ( myMesh->NbFaces() > 0 ) meshDim = 2; + if ( myMesh->NbVolumes() > 0 ) meshDim = 3; + + if ( myMeshName.empty() ) + { + int nbases = 0; + if ( cg_nbases( _fn, &nbases) == CG_OK ) + myMeshName = ( SMESH_Comment("Base_") << nbases+1 ); + else + myMeshName = "Base_0"; + } + int iBase; + if ( cg_base_write( _fn, myMeshName.c_str(), meshDim, spaceDim, &iBase ) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + if ( cg_goto(_fn, iBase, "end") != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + if ( cg_descriptor_write("About", "Created by SMESH") != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // create a structured Zone + // -------------- + + TopoDS_Shape shape = myMesh->ShapeToMesh(); + + if ( meshDim == 3 ) + { + std::set zNames; + for ( TopExp_Explorer fEx( shape, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + TopoDS_Solid currentSolid = TopoDS::Solid(fEx.Current()); + + if ( myMesh->HasStructuredGridFilled(currentSolid) ) + { + auto grid = myMesh->GetTheGrid(currentSolid).get(); + int imax = grid->nx(); + int jmax = grid->ny(); + int kmax = grid->nz(); + + cgsize_t size[9] = {imax, jmax, kmax, imax - 1, jmax - 1, kmax - 1, 0, 0, 0}; + + std::string zoneName = GetGroupName(myMesh->ShapeToIndex(currentSolid),meshDim); + zoneName = zoneName.empty() ? "ZONESOLID" + std::to_string(myMesh->ShapeToIndex(currentSolid)) : zoneName; + if ( zNames.count(zoneName) != 0 ) + zoneName = zoneName + "_" + std::to_string(myMesh->ShapeToIndex(currentSolid)); + zNames.insert( zoneName ); + // write Zone + int iZone; + if(cg_zone_write(_fn, iBase, zoneName.c_str(), size, + CGNS_ENUMV(Structured), &iZone) != CG_OK) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Grid + int iGrid=0; + if(cg_grid_write(_fn, iBase, iZone, "GridCoordinates", &iGrid) != CG_OK) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write the Coordinates + std::vector< double > coords( grid->Size() /*nx*ny*nz*/); + // write X-Coordinates + auto coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->X(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateX", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Y-Coordinates + coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->Y(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateY", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Z-Coordinates + coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->Z(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateZ", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + //end write Coordinates + + // Write Boundary condition for the grid faces + std::vector> allRanges; + std::vector boundaryNames(6,""); + grid->getAllFaceIndexLimits( allRanges ); + CheckForGroupNameOnFaceInterfaces(grid,boundaryNames); + + int faceId = 0; + std::set bNames; + for (auto pRange : allRanges) + { + if(pRange[3] < pRange[0]) {std::swap(pRange[0],pRange[3]);} + if(pRange[4] < pRange[1]) {std::swap(pRange[1],pRange[4]);} + if(pRange[5] < pRange[2]) {std::swap(pRange[2],pRange[5]);} + int cgIndexBoco = 0; + std::string boundaryName = boundaryNames[faceId].empty() ? zoneName + "_" + std::to_string(faceId) : boundaryNames[faceId]; + if ( bNames.count(boundaryName)!=0 ) + boundaryName = boundaryName + "_" + std::to_string(faceId); + bNames.insert( boundaryName ); + if(cg_boco_write(_fn, iBase, iZone, boundaryName.c_str(), CGNS_ENUMV(BCTypeNull), + CGNS_ENUMV(PointRange), 2, &pRange[0], &cgIndexBoco) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true); + faceId++; + } + // End write boundary + + // Writte Interfaces + for ( TopExp_Explorer fEx( shape, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + TopoDS_Solid neighbourSolid = TopoDS::Solid(fEx.Current()); + if ( !currentSolid.IsSame( neighbourSolid ) && myMesh->HasStructuredGridFilled(neighbourSolid)) + { + std::vector interface; + grid->GetFaceInterfaces(myMesh->GetTheGrid(neighbourSolid).get(), interface); + if ( !interface.empty() ) + { + int iConn; + std::string neigbourZoneName = GetGroupName(myMesh->ShapeToIndex(neighbourSolid),meshDim); + neigbourZoneName = neigbourZoneName.empty() ? "ZONESOLID" + std::to_string(myMesh->ShapeToIndex(neighbourSolid)) : neigbourZoneName; + std::string interfaceName = zoneName + "_" + neigbourZoneName + "_" + std::to_string(interface[0]); + if(cg_1to1_write(_fn, iBase, iZone, interfaceName.c_str(), neigbourZoneName.c_str(), + &interface[1], &interface[7], &interface[13], &iConn) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true); + } + + } + } + } + } + } + else if ( meshDim == 2 ) + { + std::set zNames; + for ( TopExp_Explorer fEx( shape, TopAbs_FACE ); fEx.More(); fEx.Next() ) + { + TopoDS_Face currentFace = TopoDS::Face(fEx.Current()); + + if ( myMesh->HasStructuredGridFilled(currentFace) ) + { + auto grid = myMesh->GetTheGrid(currentFace).get(); + int imax = grid->nx(); + int jmax = grid->ny(); + int iZone; + cgsize_t size[6] = {imax, jmax, imax - 1, jmax - 1, 0, 0}; + + std::string zoneName = GetGroupName(myMesh->ShapeToIndex(currentFace),meshDim); + zoneName = zoneName.empty() ? "ZONEFACE" + std::to_string(myMesh->ShapeToIndex(currentFace)) : zoneName; + + if ( zNames.count(zoneName) != 0 ) + zoneName = zoneName + "_" + std::to_string(myMesh->ShapeToIndex(currentFace)); + zNames.insert( zoneName ); + // write Zone + if(cg_zone_write(_fn, iBase, zoneName.c_str(), size, + CGNS_ENUMV(Structured), &iZone) != CG_OK) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Grid + int iGrid; + if(cg_grid_write(_fn, iBase, iZone, "GridCoordinates", &iGrid) != CG_OK) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write the Coordinates + std::vector< double > coords( grid->Size() /*nx*ny*nz*/); + + // write X-Coordinates + auto coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->X(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateX", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Y-Coordinates + coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->Y(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateY", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + // write Z-Coordinates + coordIter = grid->CoordinateBegin(); + for ( int i = 0; i < grid->Size(); coordIter.Next(), i++ ) + coords[i] = coordIter.Value()->Z(); + + if ( cg_coord_write( _fn, iBase, iZone, CGNS_ENUMV(RealDouble), + "CoordinateZ", &coords[0], &iGrid) != CG_OK ) + return addMessage( cg_get_error(), /*fatal = */true ); + + //End write Coordinates + + // Write Boundary condition for the grid edges + std::vector> allRanges; + std::vector boundaryNames(4,""); + grid->getAllEdgeIndexLimits( allRanges ); + CheckForGroupNameOnEdgeInterfaces(grid,boundaryNames); + int edgeId = 0; + std::set bNames; + for (auto pRange : allRanges) + { + if(pRange[2] < pRange[0]) {std::swap(pRange[0],pRange[2]);} + if(pRange[1] < pRange[3]) {std::swap(pRange[1],pRange[3]);} + int cgIndexBoco = 0; + std::string boundaryName = boundaryNames[edgeId].empty() ? zoneName + "_" + std::to_string(edgeId) : boundaryNames[edgeId]; + if ( bNames.count(boundaryName)!=0) + boundaryName = boundaryName + "_" + std::to_string(edgeId); + bNames.insert( boundaryName ); + if(cg_boco_write(_fn, iBase, iZone, boundaryName.c_str(), CGNS_ENUMV(BCTypeNull), + CGNS_ENUMV(PointRange), 2, &pRange[0], &cgIndexBoco) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true); + edgeId++; + } + // End write Boundary + + // Writte Interfaces + for ( TopExp_Explorer fEx( shape, TopAbs_FACE ); fEx.More(); fEx.Next() ) + { + TopoDS_Face neighbourFace = TopoDS::Face(fEx.Current()); + if ( !currentFace.IsSame( neighbourFace ) && myMesh->HasStructuredGridFilled(neighbourFace) ) + { + std::vector interface; + grid->GetEdgeInterfaces(myMesh->GetTheGrid(neighbourFace).get(), interface); + + if ( !interface.empty() ) + { + int iConn; + std::string neigbourZoneName = GetGroupName(myMesh->ShapeToIndex(neighbourFace),meshDim); + neigbourZoneName = neigbourZoneName.empty() ? "ZONEFACE" + std::to_string(myMesh->ShapeToIndex(neighbourFace)) : neigbourZoneName; + std::string interfaceName = zoneName + "_" + neigbourZoneName + "_" + std::to_string(interface[0]); + if(cg_1to1_write(_fn, iBase, iZone, interfaceName.c_str(), neigbourZoneName.c_str(), + &interface[1], &interface[5], &interface[9], &iConn) != CG_OK) return addMessage(cg_get_error(), /*fatal = */true); + + } + } + } + } /*end if grid filled*/ + } /*foreach face*/ + } + + return DRS_OK; +} + +//================================================================================ +/*! + * \brief Constructor + */ +//================================================================================ + +DriverStructuredCGNS_Write::DriverStructuredCGNS_Write(): _fn(0) +{ +} + +//================================================================================ +/*! + * \brief Close the cgns file at destruction + */ +//================================================================================ + +DriverStructuredCGNS_Write::~DriverStructuredCGNS_Write() +{ + if ( _fn > 0 ) + cg_close( _fn ); +} diff --git a/src/DriverCGNS/DriverStructuredCGNS_Write.hxx b/src/DriverCGNS/DriverStructuredCGNS_Write.hxx new file mode 100644 index 000000000..0198612e3 --- /dev/null +++ b/src/DriverCGNS/DriverStructuredCGNS_Write.hxx @@ -0,0 +1,70 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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 : DriverStructuredCGNS_Write.hxx +// Created : Tuesday March 19 2024 +// Author : Cesar Conopoima (cce) + +#ifndef __DriverStructuredCGNS_Write_HXX__ +#define __DriverStructuredCGNS_Write_HXX__ + +// occt +#include + +// smesh +#include "SMESH_DriverCGNS.hxx" +#include "SMESH_RegularGridTemplate.hxx" +#include "Driver_SMESHDS_Mesh.h" + +#include +#include + +/*! + * \brief Driver writinging a mesh into the CGNS file. + */ +class MESHDriverCGNS_EXPORT DriverStructuredCGNS_Write : public Driver_SMESHDS_Mesh +{ +public: + + DriverStructuredCGNS_Write(); + ~DriverStructuredCGNS_Write(); + + virtual Status Perform(); + + /*! + * \brief Search for a SMESHDS_GroupOnGeom associated to the geometry and return his name + * \param shapeToIndex, the index of the shape as given by shapeToIndex method + * \return the name of the group if found or an empty name + */ + std::string GetGroupName( const int shapeToIndex, int dim ); + + /*! + * \brief Fill name of face interfaces associated to mesh groups + */ + void CheckForGroupNameOnFaceInterfaces( const SMESHUtils::SMESH_RegularGrid* grid, std::vector& interfaceNames ); + + /*! + * \brief Fill name of edge interfaces associated to mesh groups + */ + void CheckForGroupNameOnEdgeInterfaces( const SMESHUtils::SMESH_RegularGrid* grid, std::vector& interfaceNames ); + +private: + int _fn; //!< file index +}; + +#endif diff --git a/src/SMDS/SMDS_Mesh.hxx b/src/SMDS/SMDS_Mesh.hxx index e643ee20f..50f7f5b6a 100644 --- a/src/SMDS/SMDS_Mesh.hxx +++ b/src/SMDS/SMDS_Mesh.hxx @@ -737,8 +737,7 @@ public: const int nbnodes, std::set& oldNodes ); - void setNbShapes( size_t nbShapes ); - + void setNbShapes( size_t nbShapes ); // Fields PRIVATE @@ -766,7 +765,7 @@ public: double ymin; double ymax; double zmin; - double zmax; + double zmax; }; diff --git a/src/SMESH/SMESH_Mesh.cxx b/src/SMESH/SMESH_Mesh.cxx index 2f7c6efa7..8d545cc58 100644 --- a/src/SMESH/SMESH_Mesh.cxx +++ b/src/SMESH/SMESH_Mesh.cxx @@ -53,6 +53,7 @@ #ifdef WITH_CGNS #include "DriverCGNS_Read.hxx" #include "DriverCGNS_Write.hxx" +#include "DriverStructuredCGNS_Write.hxx" #endif #include @@ -1686,6 +1687,51 @@ void SMESH_Mesh::ExportCGNS(const char * file, if ( res != Driver_Mesh::DRS_OK ) throw SALOME_Exception("Export failed"); } +//================================================================================ +/*! + * \brief Export the mesh to the StructuredCGNS file + */ +//================================================================================ + +void SMESH_Mesh::ExportStructuredCGNS(const char * file, const SMESHDS_Mesh* meshPart, const char * meshName) +{ + + int res = Driver_Mesh::DRS_OK; + SMESH_TRY; + +#ifdef WITH_CGNS + auto myMesh = meshPart ? (SMESHDS_Mesh*) meshPart : _meshDS; + + if ( myMesh->HasSomeStructuredGridFilled() ) + { + DriverStructuredCGNS_Write writer; + writer.SetFile( file ); + writer.SetMesh( const_cast( myMesh )); + writer.SetMeshName( SMESH_Comment("Mesh_") << myMesh->GetPersistentId()); + if ( meshName && meshName[0] ) + writer.SetMeshName( meshName ); + + res = writer.Perform(); + if ( res != Driver_Mesh::DRS_OK ) + { + SMESH_ComputeErrorPtr err = writer.GetError(); + if ( err && !err->IsOK() && !err->myComment.empty() ) + throw SALOME_Exception(("Export failed: " + err->myComment ).c_str() ); + } + } + +#endif + SMESH_CATCH( SMESH::throwSalomeEx ); + + if ( res == Driver_Mesh::DRS_TOO_LARGE_MESH ) + { + std::cout << "\n\n\n Going into too large mesh file path\n\n\n\n"; + throw TooLargeForExport("CGNS"); + } + + if ( res != Driver_Mesh::DRS_OK ) + throw SALOME_Exception("Export failed"); +} //================================================================================ /*! diff --git a/src/SMESH/SMESH_Mesh.hxx b/src/SMESH/SMESH_Mesh.hxx index b831998d5..f9cc2e913 100644 --- a/src/SMESH/SMESH_Mesh.hxx +++ b/src/SMESH/SMESH_Mesh.hxx @@ -305,6 +305,9 @@ class SMESH_EXPORT SMESH_Mesh const SMESHDS_Mesh* mesh, const char * meshName = 0, const bool groupElemsByType = false); + void ExportStructuredCGNS(const char * file, + const SMESHDS_Mesh* mesh, + const char * meshName = 0); void ExportGMF(const char * file, const SMESHDS_Mesh* mesh, bool withRequiredGroups = true ); diff --git a/src/SMESHDS/CMakeLists.txt b/src/SMESHDS/CMakeLists.txt index 56e87ef7b..91324502d 100644 --- a/src/SMESHDS/CMakeLists.txt +++ b/src/SMESHDS/CMakeLists.txt @@ -29,6 +29,7 @@ INCLUDE_DIRECTORIES( ${OpenCASCADE_INCLUDE_DIR} ${Boost_INCLUDE_DIRS} ${PROJECT_SOURCE_DIR}/src/SMDS + ${PROJECT_SOURCE_DIR}/src/SMESHUtils ${TBB_INCLUDES} ) @@ -47,7 +48,8 @@ SET(_link_LIBRARIES ${OpenCASCADE_KERNEL} ${OpenCASCADE_ModelingData_LIBRARIES} ${KERNEL_SALOMELocalTrace} - SMDS + SMDS + SMESHUtils ${TBB_LIBS} ) diff --git a/src/SMESHDS/SMESHDS_Mesh.cxx b/src/SMESHDS/SMESHDS_Mesh.cxx index adc79ad75..8e616d920 100644 --- a/src/SMESHDS/SMESHDS_Mesh.cxx +++ b/src/SMESHDS/SMESHDS_Mesh.cxx @@ -37,8 +37,7 @@ #include "SMESHDS_Script.hxx" #include "SMESHDS_TSubMeshHolder.hxx" -#include -#include +// #include #include #include #include @@ -48,6 +47,9 @@ #include #include +#include +#include + #include "utilities.h" class SMESHDS_Mesh::SubMeshHolder : public SMESHDS_TSubMeshHolder< const SMESHDS_SubMesh > @@ -1030,6 +1032,8 @@ void SMESHDS_Mesh::RemoveFreeElement(const SMDS_MeshElement * elt, void SMESHDS_Mesh::ClearMesh() { + myRegularGrid.Clear(); + myScript->ClearMesh(); SMDS_Mesh::Clear(); @@ -1483,6 +1487,7 @@ void SMESHDS_Mesh::SetMeshElementOnShape(const SMDS_MeshElement* anElement, //======================================================================= SMESHDS_Mesh::~SMESHDS_Mesh() { + myRegularGrid.Clear(); // myScript delete myScript; // submeshes @@ -2280,3 +2285,81 @@ bool SMESHDS_Mesh::ModifyCellNodes(vtkIdType vtkVolId, std::map localCl myGrid->ModifyCellNodes(vtkVolId, localClonedNodeIds); return true; } + +/* Create a structured grid associated to the shapeId meshed + * @remark the grid dimension (nx, ny, nz) is not associated to a space direction, those dimensions denotes which index + * run faster in the nested loop defining the nodes of a structured grid. + * for( k in range(0,nz) ) + * for( j in range(0,ny) ) + * for( i in range(0,nx) ) + * @param shape been meshed, pass the shape that would be passed to SetNodeOnFace or SetNodeInVolume method. + * @param nx, faster dimension + * @param ny, medium dimension + * @param nz, slower dimension + * @return define a new * SMESH_StructuredGrid to the shapeId. + */ +void SMESHDS_Mesh::SetStructuredGrid( const TopoDS_Shape & shape, const int nx, const int ny, const int nz ) +{ + int index = myIndexToShape.FindIndex(shape); + if ( index != 0 ) + { + if ( myRegularGrid.IsBound(index) ) + myRegularGrid.UnBind(index); + + myRegularGrid.Bind(index, std::make_shared(index,nx,ny,nz)); + } + +} + +void SMESHDS_Mesh::SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const std::shared_ptr& P, const int iIndex, const int jIndex, const int kIndex ) +{ + int index = myIndexToShape.FindIndex(shape); + if ( index != 0 && myRegularGrid.IsBound(index) ) + myRegularGrid.Seek(index)->get()->SetNode( P, iIndex, jIndex, kIndex ); +} + +void SMESHDS_Mesh::SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const SMDS_MeshNode* P, const int iIndex, const int jIndex, const int kIndex ) +{ + int index = myIndexToShape.FindIndex(shape); + if ( index != 0 && myRegularGrid.IsBound(index) ) + myRegularGrid.Seek(index)->get()->SetNode( P, iIndex, jIndex, kIndex ); +} + +void SMESHDS_Mesh::SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const SMDS_MeshNode* P, const int index ) +{ + int shapeIndex = myIndexToShape.FindIndex(shape); + if ( shapeIndex != 0 && myRegularGrid.IsBound(shapeIndex) ) + myRegularGrid.Seek(shapeIndex)->get()->SetNode( P, index ); +} + +bool SMESHDS_Mesh::HasStructuredGridFilled( const TopoDS_Shape & shape ) const +{ + int index = myIndexToShape.FindIndex(shape); + if ( index != 0 && myRegularGrid.IsBound(index) ) + return true; + else + return false; +} + +bool SMESHDS_Mesh::HasSomeStructuredGridFilled() const +{ + bool hasSomeStructuredGrid = false; + for ( TopExp_Explorer fEx( myShape, TopAbs_SOLID ); fEx.More(); fEx.Next() ) + { + TopoDS_Solid solid = TopoDS::Solid(fEx.Current()); + if ( HasStructuredGridFilled( solid ) ) hasSomeStructuredGrid = true; + } + for ( TopExp_Explorer fEx( myShape, TopAbs_FACE ); fEx.More(); fEx.Next() ) + { + TopoDS_Face face = TopoDS::Face(fEx.Current()); + if ( HasStructuredGridFilled( face ) ) hasSomeStructuredGrid = true; + } + return hasSomeStructuredGrid; +} + +const std::shared_ptr& SMESHDS_Mesh::GetTheGrid( const TopoDS_Shape & shape ) +{ + int index = myIndexToShape.FindIndex(shape); + if ( index != 0 && myRegularGrid.IsBound(index) ) + return *myRegularGrid.Seek(index); +} \ No newline at end of file diff --git a/src/SMESHDS/SMESHDS_Mesh.hxx b/src/SMESHDS/SMESHDS_Mesh.hxx index 1b8d2e2fd..9fd4c0db9 100644 --- a/src/SMESHDS/SMESHDS_Mesh.hxx +++ b/src/SMESHDS/SMESHDS_Mesh.hxx @@ -30,6 +30,8 @@ #include "SMESH_SMESHDS.hxx" #include "SMDS_Mesh.hxx" +#include "SMESH_Utils.hxx" +#include "SMESH_RegularGrid.hxx" #include "SMESHDS_SubMesh.hxx" #include @@ -679,6 +681,13 @@ class SMESHDS_EXPORT SMESHDS_Mesh : public SMDS_Mesh void CleanDownWardConnectivity(); void BuildDownWardConnectivity(bool withEdges); + virtual void SetStructuredGrid( const TopoDS_Shape & shape, const int nx, const int ny, const int nz = 1 ); + virtual void SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const std::shared_ptr& P, const int iIndex, const int jIndex, const int kIndex = 0 ); + virtual void SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const SMDS_MeshNode* point, const int iIndex, const int jIndex, const int kIndex = 0 ); + virtual void SetNodeOnStructuredGrid( const TopoDS_Shape & shape, const SMDS_MeshNode* point, const int index ); + virtual bool HasStructuredGridFilled( const TopoDS_Shape & shape ) const; + virtual bool HasSomeStructuredGridFilled() const; + virtual const std::shared_ptr& GetTheGrid( const TopoDS_Shape & shape ); ~SMESHDS_Mesh(); private: @@ -701,6 +710,9 @@ class SMESHDS_EXPORT SMESHDS_Mesh : public SMDS_Mesh int add( const SMDS_MeshElement* elem, SMESHDS_SubMesh* subMesh ); SMESHDS_SubMesh* getSubmesh( const TopoDS_Shape & shape); + + // Index the regular grid associated to the mesh in the geometry index + NCollection_DataMap> myRegularGrid; }; diff --git a/src/SMESHGUI/SMESHGUI.cxx b/src/SMESHGUI/SMESHGUI.cxx index 90f0b230d..fe6d17ac6 100644 --- a/src/SMESHGUI/SMESHGUI.cxx +++ b/src/SMESHGUI/SMESHGUI.cxx @@ -607,7 +607,7 @@ namespace { format = "CGNS"; notSupportedElemTypes.push_back( SMESH::Entity_Ball ); - } + } else if ( isGMF ) { format = "GMF"; @@ -677,6 +677,7 @@ namespace bool toFindOutDim = true; bool saveNumbers = resMgr->booleanValue( "SMESH", "med_save_numbers", true ); bool toRenumber = true; + bool structureCGNS = false; double zTol = resMgr->doubleValue( "SMESH", "med_ztolerance", 0. ); QString aFilter, aTitle = QObject::tr("SMESH_EXPORT_MESH"); @@ -697,20 +698,51 @@ namespace anInitialPath + QString("/") + aMeshName, aFilter, aTitle, false); } - else if ( isCGNS || isUNV || isDAT ) // Export to [ CGNS | UNV | DAT ] - one option + else if ( isCGNS ) + { + const char* theOptionResource = "cgns_group_elems_by_type"; + bool option = resMgr->booleanValue( "SMESH", theOptionResource, false ); + + QStringList checkBoxes; + checkBoxes << QObject::tr( "CGNS_EXPORT_ELEMS_BY_TYPE" ) << QObject::tr("STRUCTUREDCGNS"); + + SalomeApp_CheckFileDlg* fd = new SalomeApp_CheckFileDlg ( SMESHGUI::desktop(), false, checkBoxes, true, true ); + + fd->setWindowTitle( aTitle ); + fd->setNameFilter( QObject::tr( "CGNS_FILES_FILTER" ) + " (*.cgns)" ); + + if ( !anInitialPath.isEmpty() ) + fd->setDirectory( anInitialPath ); + fd->selectFile( aMeshName ); + SMESHGUI_FileValidator* fv = new SMESHGUI_FileValidator( fd ); + fd->setValidator( fv ); + fd->SetChecked( option, 0 ); + + if ( fd->exec() ) + { + aFilename = fd->selectedFile(); + structureCGNS = fd->IsChecked( 1 ); + } + + toOverwrite = fv->isOverwrite( aFilename ); + option = fd->IsChecked( 0 ); + SMESHGUI::resourceMgr()->setValue("SMESH", theOptionResource, option ); + toCreateGroups = option; + + delete fd; + } + else if ( isUNV || isDAT ) // Export to [ UNV | DAT ] - one option { - const char* theOptionResource = isCGNS ? "cgns_group_elems_by_type" : "export_renumber"; + const char* theOptionResource = "export_renumber"; bool option = resMgr->booleanValue( "SMESH", theOptionResource, false ); QStringList checkBoxes; - checkBoxes << QObject::tr( isCGNS ? "CGNS_EXPORT_ELEMS_BY_TYPE" : "SMESH_RENUMBER" ); + checkBoxes << QObject::tr( "SMESH_RENUMBER" ); SalomeApp_CheckFileDlg* fd = new SalomeApp_CheckFileDlg ( SMESHGUI::desktop(), false, checkBoxes, true, true ); fd->setWindowTitle( aTitle ); - if ( isCGNS ) - fd->setNameFilter( QObject::tr( "CGNS_FILES_FILTER" ) + " (*.cgns)" ); - else if ( isUNV ) + if ( isUNV ) fd->setNameFilter( QObject::tr( "IDEAS_FILES_FILTER" ) + " (*.unv)" ); else if ( isDAT ) fd->setNameFilter( QObject::tr( "DAT_FILES_FILTER" ) + " (*.dat)" ); @@ -726,25 +758,25 @@ namespace toOverwrite = fv->isOverwrite( aFilename ); option = fd->IsChecked( 0 ); SMESHGUI::resourceMgr()->setValue("SMESH", theOptionResource, option ); - ( isCGNS ? toCreateGroups : toRenumber ) = option; + toRenumber = option; delete fd; } else if ( isSTL ) // Export to STL { QMap aFilterMap; + QStringList filters; aFilterMap.insert( QObject::tr( "STL_ASCII_FILES_FILTER" ) + " (*.stl)", 1 ); aFilterMap.insert( QObject::tr( "STL_BIN_FILES_FILTER" ) + " (*.stl)", 0 ); - - QStringList filters; + QMap::const_iterator it = aFilterMap.begin(); for ( ; it != aFilterMap.end(); ++it ) - filters.push_back( it.key() ); - + filters.push_back( it.key() ); + SUIT_FileDlg* fd = new SUIT_FileDlg( SMESHGUI::desktop(), false, true, true ); fd->setWindowTitle( aTitle ); fd->setNameFilters( filters ); - fd->selectNameFilter( QObject::tr( "STL_ASCII_FILES_FILTER" ) + " (*.stl)" ); + if ( !anInitialPath.isEmpty() ) fd->setDirectory( anInitialPath ); fd->selectFile(aMeshName); @@ -752,7 +784,8 @@ namespace while (!is_ok) { if ( fd->exec() ) aFilename = fd->selectedFile(); - aIsASCII_STL = (aFilterMap[fd->selectedNameFilter()]) == 1 ? true: false; + if ( isSTL ) + aIsASCII_STL = (aFilterMap[fd->selectedNameFilter()]) == 1 ? true: false; is_ok = true; } delete fd; @@ -987,14 +1020,21 @@ namespace else if ( isCGNS ) { aMeshIter = aMeshList.begin(); + for( int aMeshIndex = 0; aMeshIter != aMeshList.end(); aMeshIter++, aMeshIndex++ ) { SMESH::SMESH_IDSource_var aMeshOrGroup = (*aMeshIter).first; SMESH::SMESH_Mesh_var aMeshItem = aMeshOrGroup->GetMesh(); - aMeshItem->ExportCGNS( aMeshOrGroup, - aFilename.toUtf8().data(), - toOverwrite && aMeshIndex == 0, - toCreateGroups ); + + if ( !structureCGNS ) + aMeshItem->ExportCGNS( aMeshOrGroup, + aFilename.toUtf8().data(), + toOverwrite && aMeshIndex == 0, + toCreateGroups ); + else + aMeshItem->ExportStructuredCGNS( aMeshOrGroup, + aFilename.toUtf8().data(), + toOverwrite && aMeshIndex == 0 ); } } else if ( isGMF ) @@ -1024,10 +1064,16 @@ namespace QObject::tr("SMESH_WRN_WARNING"), QObject::tr("SMESH_EXPORT_FAILED_SHORT") + "\n\n" + exText); } - else + else if ( isMED ) + { SUIT_MessageBox::warning(SMESHGUI::desktop(), QObject::tr("SMESH_WRN_WARNING"), QObject::tr("SMESH_EXPORT_FAILED") + SalomeApp_Tools::ExceptionToString(S_ex)); + } + else + SUIT_MessageBox::warning(SMESHGUI::desktop(), + QObject::tr("SMESH_WRN_WARNING"), + SalomeApp_Tools::ExceptionToString(S_ex)); wc.resume(); } } @@ -4150,7 +4196,7 @@ void SMESHGUI::createSMESHAction( const int id, const QString& po_id, const QStr pix = resMgr->loadPixmap( "SMESH", tr( QString( "ICON_%1" ).arg( po_id ).toLatin1().data() ), false ); if ( !pix.isNull() ) icon = QIcon( pix ); - + QString tooltip = tr( QString( "TOP_%1" ).arg( po_id ).toLatin1().data() ), menu = tr( QString( "MEN_%1" ).arg( po_id ).toLatin1().data() ), status_bar = tr( QString( "STB_%1" ).arg( po_id ).toLatin1().data() ); diff --git a/src/SMESHGUI/SMESH_msg_en.ts b/src/SMESHGUI/SMESH_msg_en.ts index 7d40dd113..318a047da 100644 --- a/src/SMESHGUI/SMESH_msg_en.ts +++ b/src/SMESHGUI/SMESH_msg_en.ts @@ -43,6 +43,10 @@ CGNS_EXPORT_ELEMS_BY_TYPE Group elements by type + + STRUCTUREDCGNS + Structured version (only for Hexahedron(i,j,k) and Quadrangle: Mapping meshes) + GMF_ASCII_FILES_FILTER GMF ASCII files diff --git a/src/SMESHGUI/SMESH_msg_fr.ts b/src/SMESHGUI/SMESH_msg_fr.ts index e528d3e8b..4b2af068f 100644 --- a/src/SMESHGUI/SMESH_msg_fr.ts +++ b/src/SMESHGUI/SMESH_msg_fr.ts @@ -43,6 +43,10 @@ CGNS_EXPORT_ELEMS_BY_TYPE Groupe les éléments par type + + STRUCTUREDCGNS + Format structuré (uniquement pour les maillages Hexahedron(i,j,k) et Quadrangle: Mapping) + GMF_ASCII_FILES_FILTER Fichiers GMF ASCII diff --git a/src/SMESHGUI/SMESH_msg_ja.ts b/src/SMESHGUI/SMESH_msg_ja.ts index 1cf301ad4..d7706f285 100644 --- a/src/SMESHGUI/SMESH_msg_ja.ts +++ b/src/SMESHGUI/SMESH_msg_ja.ts @@ -39,6 +39,10 @@ CGNS_EXPORT_ELEMS_BY_TYPE タイプで要素をグループ化 + + STRUCTUREDCGNS + CGNS 構造化バージョンのエクスポート + GMF_ASCII_FILES_FILTER GMFアスキーファイル diff --git a/src/SMESHUtils/CMakeLists.txt b/src/SMESHUtils/CMakeLists.txt index 4fce54a8e..872661bcc 100644 --- a/src/SMESHUtils/CMakeLists.txt +++ b/src/SMESHUtils/CMakeLists.txt @@ -63,6 +63,8 @@ SET(SMESHUtils_HEADERS SMESH_Indexer.hxx SMESH_BoostTxtArchive.hxx SMESH_MGLicenseKeyGen.hxx + SMESH_RegularGrid.hxx + SMESH_RegularGridTemplate.hxx ) # --- sources --- @@ -88,6 +90,7 @@ SET(SMESHUtils_SOURCES SMESH_PolyLine.cxx SMESH_BoostTxtArchive.cxx SMESH_MGLicenseKeyGen.cxx + SMESH_RegularGrid.cxx ) # --- rules --- diff --git a/src/SMESHUtils/SMESH_RegularGrid.cxx b/src/SMESHUtils/SMESH_RegularGrid.cxx new file mode 100644 index 000000000..60ea49bfc --- /dev/null +++ b/src/SMESHUtils/SMESH_RegularGrid.cxx @@ -0,0 +1,862 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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_StructuredGrid.cxx +// Created : Sun March 24 09:58 2024 +// Author : Cesar Conopoima (cce) + +//CAS +#include + +// SMESH +#include "SMDS_MeshNode.hxx" +#include "SMESH_TypeDefs.hxx" +#include "SMESH_RegularGridTemplate.hxx" + +using namespace SMESHUtils; + +SMESH_RegularGrid::SMESH_RegularGrid( const int id, const int nx, const int ny, const int nz ) : +myId(id), +mnx(nx), +mny(ny), +mnz(nz), +mns(nx*ny*nz) +{ + myCoordinates.Resize(0,mnx*mny*mnz-1, false); + myCoordinates.Init(nullptr); +} + +void SMESH_RegularGrid::SetNode( const std::shared_ptr& point, const int iIndex, const int jIndex, const int zIndex ) +{ + int index = zIndex * mnx * mny + jIndex * mnx + iIndex; + if ( mns > 0 && index < mns && myCoordinates[ index ] == nullptr ) + myCoordinates[ index ] = point; +} + +void SMESH_RegularGrid::SetNode( const SMDS_MeshNode* point, const int iIndex, const int jIndex, const int zIndex ) +{ + + int index = zIndex * mnx * mny + jIndex * mnx + iIndex; + if ( mns > 0 && index < mns && myCoordinates[ index ] == nullptr ) + myCoordinates[ index ] = std::move( std::make_shared(point->X(),point->Y(),point->Z())); +} + +void SMESH_RegularGrid::SetNode( const SMDS_MeshNode* point, const int index ) +{ + if ( mns > 0 && index < mns && myCoordinates[ index ] == nullptr ) + myCoordinates[ index ] = std::move( std::make_shared(point->X(),point->Y(),point->Z())); +} + +const std::shared_ptr SMESH_RegularGrid::GetNode( const int iIndex, const int jIndex, const int zIndex ) +{ + int index = zIndex * mnx * mny + jIndex * mnx + iIndex; + if ( index > -1 && index < mns ) + return myCoordinates[ index ]; + else + return nullptr; +} + + +std::tuple SMESH_RegularGrid::GetIJK(const int index) const +{ + int K = int(index/(mnx*mny)); + int J = int(index/mnx) % mny; + int I = index % mnx; + return std::tuple(I,J,K); +} + +std::vector SMESH_RegularGrid::getEdgeIndexLimitsInverted( const EdgeType edge ) const +{ + auto limits = getEdgeIndexLimits( edge ); + std::swap(limits[0],limits[2]); + std::swap(limits[1],limits[3]); + return limits; +} + +std::vector SMESH_RegularGrid::getEdgeIndexLimits( const EdgeType edge ) const +{ + switch ( edge ) { + case EdgeType::BOTTOM: + return std::vector{1,1,mnx,1}; + break; + case EdgeType::RIGHT: + return std::vector{mnx,1,mnx,mny}; + break; + case EdgeType::TOP: + return std::vector{mnx,mny,1,mny}; + break; + case EdgeType::LEFT: + return std::vector{1,mny,1,1}; + break; + default: + return std::vector{1,1,mnx,1}; + break; + } +} + +std::vector SMESH_RegularGrid::getFaceIndexLimits( const int start, const int end ) const +{ + auto startIJK = GetIJK(start); + auto endIJK = GetIJK(end); + int iStart = std::get<0>(startIJK); + int jStart = std::get<1>(startIJK); + int kStart = std::get<2>(startIJK); + int iEnd = std::get<0>(endIJK); + int jEnd = std::get<1>(endIJK); + int kEnd = std::get<2>(endIJK); + return std::vector{iStart+1,jStart+1,kStart+1,iEnd+1,jEnd+1,kEnd+1}; +} + +std::vector SMESH_RegularGrid::getEdgeIndexLimits( const int start, const int end ) const +{ + auto startIJK = GetIJK(start); + auto endIJK = GetIJK(end); + int iStart = std::get<0>(startIJK); + int jStart = std::get<1>(startIJK); + int iEnd = std::get<0>(endIJK); + int jEnd = std::get<1>(endIJK); + return std::vector{iStart+1,jStart+1,iEnd+1,jEnd+1}; +} + +std::vector SMESH_RegularGrid::getFaceIndexLimits( const FaceType face ) const +{ + switch ( face ) { + case FaceType::B_BOTTOM: + return std::vector{1,1,1,mnx,mny,1}; /*V0-V2*/ + break; + case FaceType::B_RIGHT: + return std::vector{mnx,1,1,mnx,mny,mnz}; /*V1-V6*/ + break; + case FaceType::B_BACK: + return std::vector{1,mny,1,mnx,mny,mnz}; /*V3-V6*/ + break; + case FaceType::B_LEFT: + return std::vector{1,1,1,1,mny,mnz}; /*V0-V7*/ + break; + case FaceType::B_FRONT: + return std::vector{1,1,1,mnx,1,mnz}; /*V0-V5*/ + break; + case FaceType::B_TOP: + return std::vector{1,1,mnz,mnx,mny,mnz}; /*V4-V6*/ + break; + default: + return std::vector{1,1,1,mnx,mny,1}; + break; + } +} + +int SMESH_RegularGrid::getEdgeSize( const EdgeType edge ) const +{ + switch ( edge ) { + case EdgeType::BOTTOM: + case EdgeType::TOP: + return mnx; + break; + case EdgeType::RIGHT: + case EdgeType::LEFT: + return mny; + break; + default: + return mnx; + break; + } +} + +int SMESH_RegularGrid::getFaceSize( const FaceType edge ) const +{ + switch ( edge ) { + case FaceType::B_BOTTOM: + case FaceType::B_TOP: + return mnx*mny; + break; + case FaceType::B_RIGHT: + case FaceType::B_LEFT: + return mny*mnz; + break; + case FaceType::B_FRONT: + case FaceType::B_BACK: + return mnx*mnz; + break; + default: + return mnx; + break; + } +} + + +void SMESH_RegularGrid::getAllEdgeIndexLimits(std::vector>& allRanges) +{ + this->foreachGridSide( [&]( EdgeType edge ) + { + allRanges.push_back( getEdgeIndexLimits(edge) ); + }); +} + +void SMESH_RegularGrid::getAllFaceIndexLimits(std::vector>& allRanges) +{ + this->foreachGridFace( [&]( FaceType face ) + { + allRanges.push_back( getFaceIndexLimits(face) ); + }); +} + +std::vector SMESH_RegularGrid::nodesOfFace( SMESH_RegularGrid::FaceType face ) const +{ + size_t faceSize=0; + if ( face == SMESH_RegularGrid::B_BOTTOM || face == SMESH_RegularGrid::B_TOP ) + faceSize = mnx*mny; + if ( face == SMESH_RegularGrid::B_LEFT || face == SMESH_RegularGrid::B_RIGHT ) + faceSize = mny*mnz; + if ( face == SMESH_RegularGrid::B_BACK || face == SMESH_RegularGrid::B_FRONT ) + faceSize = mnx*mnz; + + std::vector index(faceSize); + switch ( face ) { + case FaceType::B_BOTTOM: + for (size_t i = 0; i < faceSize; i++) + index[ i ] = i; + break; + case FaceType::B_TOP: + for (size_t i = 0; i < faceSize; i++) + index[ i ] = mnx * mny * mnz - mnx * mny + i; + break; + case FaceType::B_RIGHT: + for (size_t i = 0; i < faceSize; i++) + { + auto dv = std::div(i,mny); + index[ i ] = dv.quot * mnx * mny + (mnx-1)*(dv.rem+1)+dv.rem; + } + break; + case FaceType::B_LEFT: + for (size_t i = 0; i < faceSize; i++) + { + auto dv = std::div(i,mny); + index[ i ] = dv.quot * mnx * mny + (dv.rem) * mnx; + } + break; + case FaceType::B_BACK: + for (size_t i = 0; i < faceSize; i++) + { + auto dv = std::div(i,mnx); + index[ i ] = dv.quot * mnx * mny + mny * mnx - mnx + dv.rem; + } + break; + case FaceType::B_FRONT: + for (size_t i = 0; i < faceSize; i++) + { + auto dv = std::div(i,mnx); + index[ i ] = dv.quot * mnx * mny + dv.rem; + } + break; + default: + break; + } + return index; +} + +std::vector SMESH_RegularGrid::nodesOfSide( SMESH_RegularGrid::EdgeType edge ) const +{ + size_t sideSize = ( edge == SMESH_RegularGrid::BOTTOM ) || ( edge == SMESH_RegularGrid::TOP ) ? mnx : mny; + std::vector index(sideSize); + switch ( edge ) { + case EdgeType::BOTTOM: + for (size_t i = 0; i < sideSize; i++) + index[ i ] = i; + break; + case EdgeType::RIGHT: + for (size_t i = 0; i < sideSize; i++) + index[ i ] = mnx*(i+1)-1; + break; + case EdgeType::TOP: + for (size_t i = 0; i < sideSize; i++) + index[ i ] = mnx*mny-(i+1); + break; + case EdgeType::LEFT: + for (size_t i = 0; i < sideSize; i++) + index[ i ] = mnx*mny-(i+1)*mnx; + break; + default: + break; + } + return index; +} + +// Return the index of the nodes at begin and end of the edge +std::pair SMESH_RegularGrid::getEdgeLimits( const EdgeType edge ) const +{ + switch ( edge ) { + case EdgeType::BOTTOM: + return std::make_pair(0,mnx-1); + break; + case EdgeType::RIGHT: + return std::make_pair(mnx-1,mnx*mny-1); + break; + case EdgeType::TOP: + return std::make_pair(mnx*mny-1,mnx*mny-mnx); + break; + case EdgeType::LEFT: + return std::make_pair(mnx*mny-mnx,0); + break; + default: + return std::make_pair(0,mnx-1); + break; + } +} + +int SMESH_RegularGrid::getFaceCoordinateIndex( const VertexType v ) const +{ + switch ( v ) { + case VertexType::V0: + return 0; + break; + case VertexType::V1: + return mnx-1; + break; + case VertexType::V2: + return mnx*mny-1; + break; + case VertexType::V3: + return mnx*mny-mnx; + break; + case VertexType::V4: + return mnx*mny*mnz-mnx*mny; + break; + case VertexType::V5: + return mnx*mny*mnz-mnx*mny+mnx-1; + break; + case VertexType::V6: + return mnx*mny*mnz-1; + break; + case VertexType::V7: + return mnx*mny*mnz-mnx; + break; + default: + return 0; + break; + } +} + +std::tuple SMESH_RegularGrid::getFaceLimits( const FaceType face ) const +{ + switch ( face ) { + case FaceType::B_BOTTOM: /*V0-V2-V1-V3*/ + return std::tuple{getFaceCoordinateIndex(VertexType::V0), + getFaceCoordinateIndex(VertexType::V2), + getFaceCoordinateIndex(VertexType::V1), + getFaceCoordinateIndex(VertexType::V3)}; + break; + case FaceType::B_RIGHT: /*V1-V6-V2-V5*/ + return std::tuple{getFaceCoordinateIndex(VertexType::V1), + getFaceCoordinateIndex(VertexType::V6), + getFaceCoordinateIndex(VertexType::V2), + getFaceCoordinateIndex(VertexType::V5)}; + break; + case FaceType::B_BACK: /*V3-V6-V2-V7*/ + return std::tuple{getFaceCoordinateIndex(VertexType::V3), + getFaceCoordinateIndex(VertexType::V6), + getFaceCoordinateIndex(VertexType::V2), + getFaceCoordinateIndex(VertexType::V7)}; + break; + case FaceType::B_LEFT: /*V0-V7-V3-V4*/ + return std::tuple{getFaceCoordinateIndex(VertexType::V0), + getFaceCoordinateIndex(VertexType::V7), + getFaceCoordinateIndex(VertexType::V3), + getFaceCoordinateIndex(VertexType::V4)}; + break; + case FaceType::B_FRONT: /*V0-V5-V1-V4*/ + return std::tuple{getFaceCoordinateIndex(VertexType::V0), + getFaceCoordinateIndex(VertexType::V5), + getFaceCoordinateIndex(VertexType::V1), + getFaceCoordinateIndex(VertexType::V4)}; + break; + case FaceType::B_TOP: + return std::tuple{getFaceCoordinateIndex(VertexType::V4), + getFaceCoordinateIndex(VertexType::V6), + getFaceCoordinateIndex(VertexType::V5), + getFaceCoordinateIndex(VertexType::V7)}; + break; + default: + return std::tuple{getFaceCoordinateIndex(VertexType::V0), + getFaceCoordinateIndex(VertexType::V2), + getFaceCoordinateIndex(VertexType::V1), + getFaceCoordinateIndex(VertexType::V3)}; + break; + } +} + +void SMESH_RegularGrid::GetCommontInterface( FaceType face, SMESH_RegularGrid * grid, std::vector& interface ) +{ + const double tol = Precision::Confusion(); /*confusion is 1e-7, the recommended tolerance to find coincident points in 3D*/ + auto vertex = getFaceLimits( face ); + auto v0 = myCoordinates[std::get<0>(vertex)]; + auto v1 = myCoordinates[std::get<1>(vertex)]; + auto v2 = myCoordinates[std::get<2>(vertex)]; + auto v3 = myCoordinates[std::get<3>(vertex)]; + + auto foundOneTrue = [](int begin, int end, std::vector& v ) -> bool + { + for ( int i = begin; i <= end; i++ ) + if ( v[i] ) + return true; + return false; + }; + + auto foundTrueColum = [](int colum, std::vector& v ) -> bool + { + int strikeSize = 4; + for ( int i = 0; i < strikeSize; i++ ) + if ( v[colum + i*strikeSize] ) + return true; + return false; + }; + + grid->foreachGridFace( [&]( FaceType gridFace ) + { + std::vector interfaceRange,interfaceDonor,tranformationRange,tranformationDonor; + + auto neighboorVertex = grid->getFaceLimits( gridFace ); + + auto v4 = grid->GetNode(std::get<0>(neighboorVertex)); + auto v5 = grid->GetNode(std::get<1>(neighboorVertex)); + auto v6 = grid->GetNode(std::get<2>(neighboorVertex)); + auto v7 = grid->GetNode(std::get<3>(neighboorVertex)); + + + std::vector trueTable({ v0->IsEqual(*v4,tol),v0->IsEqual(*v5,tol),v0->IsEqual(*v6,tol),v0->IsEqual(*v7,tol), + v1->IsEqual(*v4,tol),v1->IsEqual(*v5,tol),v1->IsEqual(*v6,tol),v1->IsEqual(*v7,tol), + v2->IsEqual(*v4,tol),v2->IsEqual(*v5,tol),v2->IsEqual(*v6,tol),v2->IsEqual(*v7,tol), + v3->IsEqual(*v4,tol),v3->IsEqual(*v5,tol),v3->IsEqual(*v6,tol),v3->IsEqual(*v7,tol)}); + + std::vector trueCorner({foundOneTrue(0,3,trueTable), + foundOneTrue(4,7,trueTable), + foundOneTrue(8,11,trueTable), + foundOneTrue(12,15,trueTable) }); + + if ( trueCorner[0] && + trueCorner[1] && + trueCorner[2] && + trueCorner[3] ) /*Face to Face interface 100% conform*/ + { + interfaceRange = this->getFaceIndexLimits( face ); + interfaceDonor = grid->getFaceIndexLimits( gridFace ); + } + else if ( trueCorner[0] || trueCorner[1] || + trueCorner[2] || trueCorner[3] ) /*Partial Face to Face. Only one intersection then all the other 3 vertex are: 2 in the edges 1 inside the face*/ + { + // Two possible cases: + // 1) face < gridFace + if ( this->getFaceSize( face ) < grid->getFaceSize( gridFace ) ) + { + auto nodeToSearch = trueCorner[0] ? v1 : trueCorner[2] ? v3 : trueCorner[1] ? v0 : v2; + + grid->foreachNodeOnFace( gridFace, [&] (const std::shared_ptr sidePoint, const int nodeIndex) + { + if ( nodeToSearch->IsEqual( *sidePoint, tol ) ) + { + interfaceRange = this->getFaceIndexLimits( face ); + + auto startIndex = foundTrueColum(0,trueTable) ? std::get<0>(neighboorVertex) : + foundTrueColum(2,trueTable) ? std::get<2>(neighboorVertex) : + foundTrueColum(1,trueTable) ? std::get<1>(neighboorVertex) : + std::get<3>(neighboorVertex); + + interfaceDonor = startIndex < nodeIndex ? + grid->getFaceIndexLimits( startIndex, nodeIndex ) : + grid->getFaceIndexLimits( nodeIndex, startIndex ); + return; + } + }); + } + // 2) face > gridEdge + else if ( this->getFaceSize( face ) > grid->getFaceSize( gridFace ) ) + { + auto nodeToSearch = foundTrueColum(0,trueTable) ? v5 : + foundTrueColum(2,trueTable) ? v7 : + foundTrueColum(1,trueTable) ? v4 : v6; + + this->foreachNodeOnFace( face, [&] (const std::shared_ptr sidePoint, const int nodeIndex ) + { + if ( nodeToSearch->IsEqual( *sidePoint, tol ) ) + { + auto startIndex = trueCorner[0] ? std::get<0>(vertex) : + trueCorner[2] ? std::get<2>(vertex) : + trueCorner[1] ? std::get<1>(vertex) : std::get<3>(vertex); + + interfaceRange = startIndex < nodeIndex ? + this->getFaceIndexLimits( startIndex, nodeIndex) : + this->getFaceIndexLimits( nodeIndex, startIndex ); + + interfaceDonor = grid->getFaceIndexLimits( gridFace ); + return; + } + }); + } + } + + if ( !interfaceRange.empty() && !interfaceDonor.empty() ) + { + tranformationRange = this->computeTransformation( face, gridFace, interfaceRange, interfaceDonor ); + tranformationDonor = grid->computeTransformation( gridFace, face, interfaceDonor, interfaceRange ); + interface = std::vector{(int)face}+interfaceRange+interfaceDonor+tranformationRange; + auto dualInterface = std::vector{(int)face}+interfaceDonor+interfaceRange+tranformationDonor; + this->setInterface( face, grid->id(), interface ); + grid->setInterface( gridFace, this->id(), dualInterface ); + } + }); +} + +void SMESH_RegularGrid::GetCommontInterface( EdgeType edge, SMESH_RegularGrid * grid, std::vector& interface ) +{ + const double tol = Precision::Confusion(); /*confusion is 1e-7, the recommended tolerance to find coincident points in 3D*/ + auto vertex = getEdgeLimits( edge ); + auto v0 = myCoordinates[vertex.first]; + auto v1 = myCoordinates[vertex.second]; + + grid->foreachGridSide( [&]( EdgeType gridEdge ) + { + std::vector interfaceRange,interfaceDonor,tranformationRange,tranformationDonor; + auto neighboorVertex = grid->getEdgeLimits( gridEdge ); + auto v2 = grid->GetNode(neighboorVertex.first); + auto v3 = grid->GetNode(neighboorVertex.second); + + std::vector trueTable({ v0->IsEqual(*v2,tol), + v0->IsEqual(*v3,tol), + v1->IsEqual(*v2,tol), + v1->IsEqual(*v3,tol) }); + + if ( (trueTable[0] || trueTable[1] ) && /*Case start-end vertex are coincident in both edges trivial*/ + (trueTable[2] ||trueTable[3] ) ) + { + interfaceRange = this->getEdgeIndexLimits( edge ); + interfaceDonor = grid->getEdgeIndexLimitsInverted( gridEdge ); + } + else if ( trueTable[0] || /*Case start OR end vertex are coincident in both edges*/ + trueTable[1] || + trueTable[2] || + trueTable[3] ) + { + // Two possible cases: + // 1) edge < gridEdge means that v1 is going to be found as an intersection inside gridEdge side!!! + if ( this->getEdgeSize( edge ) < grid->getEdgeSize( gridEdge ) ) + { + auto nodeToSearch = (trueTable[0] || trueTable[1]) ? v1 : v0; + grid->foreachNodeOnSide( gridEdge, [&] (const std::shared_ptr sidePoint, const int nodeIndex) + { + if ( nodeToSearch->IsEqual( *sidePoint, tol ) ) + { + interfaceRange = this->getEdgeIndexLimits( edge ); + auto startIndex = (trueTable[0] || trueTable[2]) ? neighboorVertex.first : neighboorVertex.second; + interfaceDonor = startIndex < nodeIndex ? + grid->getEdgeIndexLimits( startIndex, nodeIndex ) : + grid->getEdgeIndexLimits( nodeIndex, startIndex ); + + } + }); + } + // 2) edge > gridEdge means that v2 or v3 are going to be found as an intersection inside edge side!!! + else if ( this->getEdgeSize( edge ) > grid->getEdgeSize( gridEdge ) ) + { + auto nodeToSearch = (trueTable[0] || trueTable[2]) ? v3 : v2; + this->foreachNodeOnSide( edge, [&] (const std::shared_ptr sidePoint, const int nodeIndex) + { + if ( nodeToSearch->IsEqual( *sidePoint, tol ) ) + { + auto startIndex = (trueTable[0] || trueTable[1]) ? vertex.first : vertex.second; + interfaceRange = startIndex < nodeIndex ? + this->getEdgeIndexLimits( startIndex, nodeIndex) : + this->getEdgeIndexLimits( nodeIndex, startIndex ); + + interfaceDonor = grid->getEdgeIndexLimitsInverted( gridEdge ); + } + }); + } + } + else /*Case nor start nor end are coincident in both edges but one of the edges can be embemded completly in the other*/ + { + // Discard quickly non parallel and far edges + auto v01 = v1->Coord()-v0->Coord(); + auto v23 = v3->Coord()-v2->Coord(); + if ( std::fabs( v01.Dot(v23) ) == 1.0 /*parallel edges might have intersection*/ ) + { + double v01Module = v01.Modulus(); + double v23Module = v23.Modulus(); + auto ref = v01Module > v23Module ? v01 : v23; + auto v0Ref = v01Module > v23Module ? v0->Coord() : v2->Coord(); + auto v2Ref = v01Module > v23Module ? std::pair(v2->Coord(),v3->Coord()) : std::pair(v0->Coord(),v1->Coord()); + + auto t = ref.X() != 0. ? std::pair( (v2Ref.first.X()-v0Ref.X())/ref.X(), (v2Ref.second.X()-v0Ref.X())/ref.X()) : + ref.Y() != 0. ? std::pair( (v2Ref.first.Y()-v0Ref.Y())/ref.Y(), (v2Ref.second.Y()-v0Ref.Y())/ref.Y()) : + std::pair( (v2Ref.first.Z()-v0Ref.Z())/ref.Z(), (v2Ref.second.Z()-v0Ref.Z())/ref.Z()); + + // check if the smaller edge contains a point inside the bigger one + // use the 3D line equation for this. + // 1) t goes from [0,1] in the reference (bigger) edge. Check if any of the points of the other edge is 0 <= t12 <= 1 + + if ( (t.first >= 0. && t.first <= 1.) || (t.second >= 0. && t.second <= 1.) ) + { + bool isFirstVertex = t.first >= 0. && t.first <= 1.; + bool isSecondVertex = t.second >= 0. && t.second <= 1.; + + if ( v01Module > v23Module && (isFirstVertex && isSecondVertex) /*100% edge embeded in the other && edge>gridEdge*/ ) + { + int startIndex = -1; + int endIndex = -1; + bool firstToBeFound = true; + this->foreachNodeOnSide( edge, [&] (const std::shared_ptr sidePoint, const int nodeIndex) + { + if ( v2->IsEqual( *sidePoint, tol ) || v3->IsEqual( *sidePoint, tol ) ) + { + if ( firstToBeFound ) + { + startIndex = nodeIndex; + firstToBeFound = false; + } + else + endIndex = nodeIndex; + } + }); + + if ( startIndex != -1 && endIndex != -1 ) + { + interfaceRange = this->getEdgeIndexLimits( startIndex, endIndex ); + interfaceDonor = grid->getEdgeIndexLimitsInverted( gridEdge ); + } + } + else if ( (isFirstVertex && isSecondVertex) /*100% edge embeded in the other && gridEdge>edge*/ ) + { + int startIndex = -1; + int endIndex = -1; + bool firstToBeFound = true; + grid->foreachNodeOnSide( gridEdge, [&] (const std::shared_ptr sidePoint, const int nodeIndex) + { + if ( v0->IsEqual( *sidePoint, tol ) || v1->IsEqual( *sidePoint, tol ) ) + { + if ( firstToBeFound ) + { + startIndex = nodeIndex; + firstToBeFound = false; + } + else + endIndex = nodeIndex; + } + }); + + if ( startIndex != -1 && endIndex != -1 ) + { + interfaceRange = this->getEdgeIndexLimits( edge ); + interfaceDonor = grid->getEdgeIndexLimits( startIndex, endIndex ); + } + } + } + } + } + + if ( !interfaceRange.empty() && !interfaceDonor.empty() ) + { + tranformationRange = this->computeTransformation( edge, gridEdge, interfaceRange, interfaceDonor ); + tranformationDonor = grid->computeTransformation( gridEdge, edge, interfaceDonor, interfaceRange ); + interface = std::vector{(int)edge} + interfaceRange+interfaceDonor+tranformationRange; + auto dualInterface = std::vector{(int)edge} + interfaceDonor+interfaceRange+tranformationDonor; + this->setInterface( edge, grid->id(), interface ); + grid->setInterface( gridEdge, this->id(), dualInterface ); + } + }); +} + +std::vector SMESH_RegularGrid::computeTransformation( const SMESH_RegularGrid::EdgeType edge, SMESH_RegularGrid::EdgeType gridDonorEdge, + std::vector& interfaceRange, std::vector& interfaceDonor ) const +{ + std::vector transform(2); + int r[2], d[2]; + for(int i = 0; i < 2; i++) { + r[i] = interfaceRange[i + 2] - interfaceRange[i]; + d[i] = interfaceDonor[i + 2] - interfaceDonor[i]; + } + for(int i = 0; i < 2; i++) { + transform[i] = 0; + for(int j = 0; j < 2; j++) { + if(std::abs(r[i]) == std::abs(d[j])) { // == 0 on an interface + transform[i] = j + 1; + if(!r[i] && !d[j]) { // on an interface + // both interfaces correspond to a min index or to a max index + bool topBottom = (edge == SMESH_RegularGrid::BOTTOM || edge == SMESH_RegularGrid::TOP) && (gridDonorEdge == SMESH_RegularGrid::BOTTOM || gridDonorEdge == SMESH_RegularGrid::TOP); + bool leftRight = (edge == SMESH_RegularGrid::RIGHT || edge == SMESH_RegularGrid::LEFT) && (gridDonorEdge == SMESH_RegularGrid::RIGHT || gridDonorEdge == SMESH_RegularGrid::LEFT); + if( topBottom || leftRight ) + transform[i] *= -1; + } + else { + if(r[i] * d[j] < 0) transform[i] *= -1; + } + } + } + } + return transform; +} + +std::vector SMESH_RegularGrid::computeTransformation( const SMESH_RegularGrid::FaceType face, SMESH_RegularGrid::FaceType gridDonorFace, + std::vector& interfaceRange, std::vector& interfaceDonor ) const +{ + std::vector transform(3); + int r[3], d[3]; + for(int i = 0; i < 3; i++) { + r[i] = interfaceRange[i + 3] - interfaceRange[i]; + d[i] = interfaceDonor[i + 3] - interfaceDonor[i]; + } + for(int i = 0; i < 3; i++) { + transform[i] = 0; + for(int j = 0; j < 3; j++) { + if(std::abs(r[i]) == std::abs(d[j])) { // == 0 on an interface + transform[i] = j + 1; + if(!r[i] && !d[j]) { // on an interface + // both interfaces correspond to a min index or to a max index + bool bothMin = ( face == SMESH_RegularGrid::B_BOTTOM || + face == SMESH_RegularGrid::B_LEFT || + face == SMESH_RegularGrid::B_FRONT ) && + (gridDonorFace == SMESH_RegularGrid::B_BOTTOM || + gridDonorFace == SMESH_RegularGrid::B_LEFT || + gridDonorFace == SMESH_RegularGrid::B_FRONT ); + + bool bothMax = (face == SMESH_RegularGrid::B_TOP || + face == SMESH_RegularGrid::B_RIGHT || + face == SMESH_RegularGrid::B_BACK ) && + (gridDonorFace == SMESH_RegularGrid::B_TOP || + gridDonorFace == SMESH_RegularGrid::B_RIGHT || + gridDonorFace == SMESH_RegularGrid::B_FRONT ); + + if( bothMin || bothMax ) + transform[i] *= -1; + } + else { + if(r[i] * d[j] < 0) transform[i] *= -1; + } + } + } + } + return transform; +} + +void SMESH_RegularGrid::setInterface( SMESH_RegularGrid::EdgeType edge, const int gridId, std::vector& interface ) +{ + myInterfaceMap[edge][gridId] = interface; +} + +void SMESH_RegularGrid::setInterface( SMESH_RegularGrid::FaceType face, const int gridId, std::vector& interface ) +{ + myFaceInterafaceMap[face][gridId] = interface; +} + +SMESH_RegularGrid::FaceType SMESH_RegularGrid::getFaceTypeByGeomFace( TopoDS_Shape shapeFace ) const +{ + const double tol = Precision::Confusion(); + SMESH_RegularGrid::FaceType foundFace = SMESH_RegularGrid::B_NONE; + this->foreachGridFace( [&]( FaceType face ) + { + auto faceLimits = getFaceLimits(face); + auto v0 = myCoordinates[std::get<0>(faceLimits)]; + auto v1 = myCoordinates[std::get<1>(faceLimits)]; + auto v2 = myCoordinates[std::get<2>(faceLimits)]; + auto v3 = myCoordinates[std::get<3>(faceLimits)]; + int nVertex = 0; + int nVertexInFace = 0; + for ( TopExp_Explorer fEx( shapeFace, TopAbs_VERTEX ); fEx.More(); fEx.Next() ) + { + gp_Pnt pV = BRep_Tool::Pnt( TopoDS::Vertex(fEx.Current()) ); + nVertexInFace += int( v0->IsEqual(pV,tol) ) + int( v1->IsEqual(pV,tol) ) + int( v2->IsEqual(pV,tol) ) + int( v3->IsEqual(pV,tol) ) ; + nVertex++; + } + if ( nVertex == nVertexInFace ) + { + foundFace = face; + return; + } + }); + + return foundFace; +} + +SMESH_RegularGrid::EdgeType SMESH_RegularGrid::getEdgeTypeByGeomEdge( TopoDS_Shape shapeEdge ) const +{ + const double tol = Precision::Confusion(); + SMESH_RegularGrid::EdgeType foundEdge = SMESH_RegularGrid::NONE; + this->foreachGridSide ( [&]( EdgeType edge ) + { + auto edgeLimits = getEdgeLimits(edge); + auto v0 = myCoordinates[edgeLimits.first]; + auto v1 = myCoordinates[edgeLimits.second]; + int nVertex = 0; + int nVertexInEdge = 0; + for ( TopExp_Explorer fEx( shapeEdge, TopAbs_VERTEX ); fEx.More(); fEx.Next() ) + { + gp_Pnt pV = BRep_Tool::Pnt( TopoDS::Vertex(fEx.Current()) ); + nVertexInEdge += int( v0->IsEqual(pV,tol) ) + int( v1->IsEqual(pV,tol) ); + nVertex++; + } + if ( nVertex == nVertexInEdge ) + { + foundEdge = edge; + return; + } + }); + + return foundEdge; +} + +bool SMESH_RegularGrid::GetPrecomputedInterface( EdgeType edge, const SMESH_RegularGrid * grid, std::vector& interface ) +{ + if ( myInterfaceMap.count(edge) && myInterfaceMap[edge].count(grid->id()) ) + { + interface = myInterfaceMap[edge][grid->id()]; + return true; + } + return false; +} + +bool SMESH_RegularGrid::GetPrecomputedInterface( FaceType face, const SMESH_RegularGrid * grid, std::vector& interface ) +{ + if ( myFaceInterafaceMap.count(face) && myFaceInterafaceMap[face].count(grid->id()) ) + { + interface = myFaceInterafaceMap[face][grid->id()]; + return true; + } + return false; +} + +void SMESH_RegularGrid::GetEdgeInterfaces( SMESH_RegularGrid * grid, std::vector& interface ) +{ + this->foreachGridSide( [&]( EdgeType edge ) + { + if ( !this->GetPrecomputedInterface( edge, grid, interface ) ) + this->GetCommontInterface( edge, grid, interface ); + }); +} + +void SMESH_RegularGrid::GetFaceInterfaces( SMESH_RegularGrid * grid, std::vector& interface ) +{ + this->foreachGridFace( [&]( FaceType face ) + { + if ( !this->GetPrecomputedInterface( face, grid, interface ) ) + this->GetCommontInterface( face, grid, interface ); + }); +} + + +SMESH_RegularGrid::~SMESH_RegularGrid() +{ + myInterfaceMap.clear(); + myFaceInterafaceMap.clear(); +} \ No newline at end of file diff --git a/src/SMESHUtils/SMESH_RegularGrid.hxx b/src/SMESHUtils/SMESH_RegularGrid.hxx new file mode 100644 index 000000000..2a437ecf6 --- /dev/null +++ b/src/SMESHUtils/SMESH_RegularGrid.hxx @@ -0,0 +1,363 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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_RegularGrid.hxx +// Created : Sun March 24 09:58 2024 +// Author : Cesar Conopoima (cce) + + +#ifndef __SMESH_RegularGrid_HXX__ +#define __SMESH_RegularGrid_HXX__ + +//OCC +#include +#include + +//STD +#include +#include + +#include "SMESH_Utils.hxx" +#include "SMDS_MeshNode.hxx" + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +class Adaptor3d_Surface; +class Adaptor2d_Curve2d; +class Adaptor3d_Curve; +class gp_Pnt; +class gp_XYZ; + +template +std::vector operator+(std::vector const &x, std::vector const &y) +{ + std::vector vec; + vec.reserve(x.size() + y.size()); + vec.insert(vec.end(), x.begin(), x.end()); + vec.insert(vec.end(), y.begin(), y.end()); + return vec; +} + +namespace SMESHUtils +{ + /*! + * \brief Define a regular grid of nx,ny,nz dimension. A vector with gp_Pnt is + * defined by default by a canonical order, the canonical order is: + * i-index goes faster, then j-index and finaly z-index. + * This data structure was created to index mesh nodes created by Quadrangle_2D and Prism_3D meshers. + * The StructuredCGNS export driver uses this data structure to define zones, grid coordinate points and interfaces. + */ + class SMESHUtils_EXPORT SMESH_RegularGrid + { + public: + + SMESH_RegularGrid( const int id, const int nx, const int ny, const int nz = 1 ); + void SetNode( const std::shared_ptr& point, const int iIndex, const int jIndex, const int zIndex = 0 ); + void SetNode( const SMDS_MeshNode* point, const int iIndex, const int jIndex, const int zIndex = 0 ); + void SetNode( const SMDS_MeshNode* point, const int index ); + const std::shared_ptr GetNode( const int index ) const { return myCoordinates[index]; }; + const std::shared_ptr GetNode( const int iIndex, const int jIndex, const int zIndex ); + int nx() const { return mnx; }; + int ny() const { return mny; }; + int nz() const { return mnz; }; + int Size() { return mns;}; + NCollection_Array1>::Iterator CoordinateBegin() { return NCollection_Array1>::Iterator( myCoordinates ); }; + int id() const { return myId; }; + + /* + * + * o------TOP (j=ny-1)----o + * | | + * | | + * | | + * LEFT (i=0) RIGTH(i=nx-1) + * | | + * | | + * | | + * o------BOTTOM(j=0)-----o + */ + enum EdgeType + { + BOTTOM, + RIGHT, + TOP, + LEFT, + NONE /**/ + }; + + enum VertexType + { + V0, + V1, + V2, + V3, + V4, + V5, + V6, + V7 + }; + /* Define vertex numeration convention + * + * V7 V6 + * +--------+--------+--------+ + * / / / /| + * +--------+--------+--------+ | + * / / / /| | + * +--------+--------+--------+ | + + * / / / V5 /| |/| + * V4 +--------+--------+--------+ | + | + * | | | |/| | + * | | | + | + + * | | |/| |/| + * + | | + | + * | | | |/| | + * | +-------+----------+-| + | + V2 + * | / V3 |/| |/ + * + / | + + * | / | |/ + * | / | + + * |/ |/ + * +--------+--------+--------+ + * V0 V1 + * + * Canonical cartesian axis orientation + * + * ^ ^ j (or y) + * k (or z)| / + * | / + * |/ + * +------> i (or x) + * + * + */ + + enum FaceType + { + B_BOTTOM, /*k=0 move (i,j)*/ + B_RIGHT, /*i=mnx move (j,k)*/ + B_BACK, /*j=mny move (i,k)*/ + B_LEFT, /*i=0 move (j,k) */ + B_FRONT, /*j=0 move (i,k)*/ + B_TOP, /*k=mnz move (i,j) */ + B_NONE /**/ + }; + + /* + * +--------+--------+--------+ + * / / / /| + * +--------+--------+--------+ | + * / / Top / /| | + * +--------+--------+--------+ | + + * / / / /| |/| + * +--------+--------+--------+ | + | + * | | | | |/| | <- Right + * | | | | + | + + * | | | |/| |/| + * +--------+--------+--------+ | + | + * Left -> | | Front | | |/| | + * | | | | + | + + * | | | |/| |/ + * +--------+--------+--------+ | + + * | | | | |/ + * | | | | + + * | | | |/ + * +--------+--------+--------+ + * ^ + * | + * Bottom + * + * Canonical cartesian axis orientation + * + * ^ ^ j (or y) + * k (or z)| / + * | / + * |/ + * +------> i (or x) + * + * + */ + + // \brief Given a grid and a potentially neigbor grid return the vector describing the interface used by the cgns mesh format + // \remark the interface vector has the follow informations: + // id the target range the donor range the transformation + // {edgeId, ibegin,jbegin, iend,jend, transformation} + // \remark As interfaces are unique, the dual interface (interaction of the neighbor grid and the current one) is also defined + // \param the neighbor grid + // \return the interface vector filled + void GetEdgeInterfaces( SMESH_RegularGrid * grid, std::vector& interface ); + + // \brief Given a grid and a potentially neigbor grid return the vector describing the interface used by the cgns mesh format + // \remark the interface vector has the follow informations: + // id the target range the donor range the transformation + // {edgeId, ibegin,jbegin,kbegin, iend,jend,kend, transformation} + // \remark As interfaces are unique, the dual interface (interaction of the neighbor grid and the current one) is also defined + // \param the neighbor grid + // \return the interface vector filled + void GetFaceInterfaces( SMESH_RegularGrid * grid, std::vector& interface ); + + // \brief Fill the allRanges vector with the boundaries of the grid + void getAllEdgeIndexLimits( std::vector>& allRanges ); + + // \brief Fill the allRanges vector with the boundaries of the grid + void getAllFaceIndexLimits( std::vector>& allRanges ); + + // \brief Get limits of the edge in the order (ibegin,jbegin,iend,jend)+1 because index in CGNS are not zero based + std::vector getEdgeIndexLimits( const EdgeType edge ) const; + + // \brief Get limits of the edge in the order (iend,jend,ibegin,jbegin)+1 because index in CGNS are not zero based + std::vector getEdgeIndexLimitsInverted( const EdgeType edge ) const; + + // \brief Get limits of the face in the order (ibegin,jbegin,kbegin,iend,jend,kend)+1 because index in CGNS are not zero based + std::vector getFaceIndexLimits( SMESH_RegularGrid::FaceType face ) const; + + // \brief Return the faceType to which the passed geometrical face belows to + SMESH_RegularGrid::FaceType getFaceTypeByGeomFace( TopoDS_Shape shapeFace ) const; + + // \brief Return the edgeType to which the passed geometrical edge belows to + SMESH_RegularGrid::EdgeType getEdgeTypeByGeomEdge( TopoDS_Shape shapeEdge ) const; + + ~SMESH_RegularGrid(); + + protected: + + // Utility functions + + // \brief Given one of the edge types and a grid determine whether or not they are adjacent on that side + // \remark There are three types of interfaces. + // 1) end to end intersection, both sides are exactly the same + // 2) one of the sides is smaller than the other but the grid are conform so the intersection is given exacly node by node + // 3) the same for 2) but nodes are not coincident + // \param edge, the edge where we are checking for interface + // \param grid, a neighbor grid + // \return interface, the interface vector filled in case an interface was found + void GetCommontInterface( EdgeType edge, SMESH_RegularGrid * grid, std::vector& interface ); + void GetCommontInterface( FaceType face, SMESH_RegularGrid * grid, std::vector& interface ); + + // \brief Accessor method to interfaces already computed by other grid + // \param edge, the edge where we are checking for interface + // \param grid, a neighbor grid + // \param interface, the interface vector filled in case an interface was found + bool GetPrecomputedInterface( EdgeType edge, const SMESH_RegularGrid * grid, std::vector& interface ); + + // \brief Accessor method to interfaces already computed by other grid + // \param face, the face where we are checking for interface + // \param grid, a neighbor grid + // \param interface, the interface vector filled in case an interface was found + bool GetPrecomputedInterface( FaceType face, const SMESH_RegularGrid * grid, std::vector& interface ); + + // \brief Set the interface found between two neighbor grids. + // \remark fill the inner map between edges and grids used to recover precomputed interfaces with GetPrecomputedInterface method + void setInterface( SMESH_RegularGrid::EdgeType edge, const int gridId, std::vector& interface ); + + // \brief Set the interface found between two neighbor grids. + // \remark fill the inner map between edges and grids used to recover precomputed interfaces with GetPrecomputedInterface method + void setInterface( SMESH_RegularGrid::FaceType face, const int gridId, std::vector& interface ); + + // \brief Get the coordinate indexes defining the limit of the edge + std::pair getEdgeLimits( const SMESH_RegularGrid::EdgeType edge ) const; + + // \brief Get the coordinate indexes defining the limit of the face + std::tuple getFaceLimits( const SMESH_RegularGrid::FaceType face ) const; + + // \brief Get limits of the edge in the order (ibegin,jbegin,iend,jend)+1 because index in CGNS are not zero based + std::vector getEdgeIndexLimits( const int start, const int end ) const; + std::vector getFaceIndexLimits( const int start, const int end ) const; + + // \brief Compute the tranformation vector following the rules of the cgns format + std::vector computeTransformation( const SMESH_RegularGrid::EdgeType edge, SMESH_RegularGrid::EdgeType gridDonorEdge, std::vector& interfaceRange, std::vector& interfaceDonor ) const; + + // \brief Compute the tranformation vector following the rules of the cgns format + std::vector computeTransformation( const SMESH_RegularGrid::FaceType face, SMESH_RegularGrid::FaceType gridDonorFace, std::vector& interfaceRange, std::vector& interfaceDonor ) const; + + // \brief Number of nodes on the edge + int getEdgeSize(const EdgeType edge ) const; + + // \brief Number of nodes on the face + int getFaceSize(const FaceType face ) const; + + // \brief Coordinate index of the given vertex + int getFaceCoordinateIndex( const VertexType v ) const; + + // \brief Given the index of a coordinate return his (i,j,k) location in the grid + // \param the index in the coordinates vector + // \return a tuple with the (i,j,k) coordinates of the point in the grid + std::tuple GetIJK( const int index ) const; + + // \brief Iterator function to execute the 'fSide' function for each edge side of the grid + // \remark utility function to be used as: + // grid->foreachGridSide( [&]( EdgeType edge ) + // { + // // Inner logic to be execute for each side of the grid + // }); + template + void foreachGridSide( const FUNCTION& fSide ) const; + + // \brief Iterator function to execute the 'fSide' function for each face side of the grid + // \remark utility function to be used as: + // grid->foreachGridFace( [&]( FaceType face ) + // { + // // Inner logic to be execute for each face side of the grid + // }); + template + void foreachGridFace( const FUNCTION& fSide ) const; + + // \brief Iterator function to execute the lambda function with argument for each grid point and index in the gridEdge + // \remark utility function to be used as: + // grid->foreachNodeOnSide( gridEdge, [&]( (const std::shared_ptr , const int nodeIndex ) ) + // { + // // Inner logic to be execute for each node and index of the passed grid edge side + // }); + void foreachNodeOnSide( SMESH_RegularGrid::EdgeType edge, const std::function point, const int index )>& function ) const; + + // \brief Iterator function to execute the lambda function with argument for each grid point and index in the gridFace + // \remark utility function to be used as: + // grid->foreachNodeOnFace( gridFace, [&]( (const std::shared_ptr, const int nodeIndex ) ) + // { + // // Inner logic to be execute for each node and index of the passed grid face side + // }); + void foreachNodeOnFace( SMESH_RegularGrid::FaceType face, const std::function point, const int index )>& function ) const; + + // \brief Get all the indices of nodes in an edge side + // \param edge, the edge of the grid + // \return the index of all the nodes in the edge + std::vector nodesOfSide( SMESH_RegularGrid::EdgeType edge ) const; + + // \brief Get all the indices of nodes in an face side + // \param face, the face of the grid + // \return the index of all the nodes in the face + std::vector nodesOfFace( SMESH_RegularGrid::FaceType face ) const; + + private: + int myId; + int mnx,mny,mnz,mns; + NCollection_Array1> myCoordinates; + + std::map>> myInterfaceMap; + std::map>> myFaceInterafaceMap; + }; +} + +#endif diff --git a/src/SMESHUtils/SMESH_RegularGridTemplate.hxx b/src/SMESHUtils/SMESH_RegularGridTemplate.hxx new file mode 100644 index 000000000..64972decc --- /dev/null +++ b/src/SMESHUtils/SMESH_RegularGridTemplate.hxx @@ -0,0 +1,63 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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_TemplateRegularGrid.hxx +// Created : Fry April 05 11:20 2024 +// Author : Cesar Conopoima (cce) + + +#ifndef __SMESH_RegularGridTemplate_HXX__ +#define __SMESH_RegularGridTemplate_HXX__ + +//STD +#include + +//SMESH::Utils +#include "SMESH_RegularGrid.hxx" + +namespace SMESHUtils +{ + template + void SMESH_RegularGrid::foreachGridSide( const FUNCTION& fSide ) const + { + for (size_t edgeType = BOTTOM; edgeType <= LEFT; edgeType++) + fSide( static_cast( edgeType ) ); + } + + template + void SMESH_RegularGrid::foreachGridFace( const FUNCTION& fSide ) const + { + for (size_t faceType = B_BOTTOM; faceType <= B_TOP; faceType++) + fSide( static_cast( faceType ) ); + } + + void SMESH_RegularGrid::foreachNodeOnSide( SMESH_RegularGrid::EdgeType edge, const std::function point, const int index )>& function ) const + { + for (auto node : nodesOfSide(edge)) + function(myCoordinates[node],node); + } + + void SMESH_RegularGrid::foreachNodeOnFace( SMESH_RegularGrid::FaceType face, const std::function point, const int index )>& function ) const + { + for (auto node : nodesOfFace(face)) + function(myCoordinates[node],node); + } +} + + +#endif diff --git a/src/SMESH_I/SMESH_Mesh_i.cxx b/src/SMESH_I/SMESH_Mesh_i.cxx index e1c8ca916..7f732ae9c 100644 --- a/src/SMESH_I/SMESH_Mesh_i.cxx +++ b/src/SMESH_I/SMESH_Mesh_i.cxx @@ -4511,7 +4511,7 @@ void SMESH_Mesh_i::ExportPartToSTL(::SMESH::SMESH_IDSource_ptr meshPart, //================================================================================ /*! - * \brief Export a part of mesh to an STL file + * \brief Export a part of mesh to an CGNS file */ //================================================================================ @@ -4551,6 +4551,50 @@ void SMESH_Mesh_i::ExportCGNS(::SMESH::SMESH_IDSource_ptr meshPart, #endif } +//================================================================================ +/*! + * \brief Export a part of mesh to an StructuredCGNS file + */ +//================================================================================ + +void SMESH_Mesh_i::ExportStructuredCGNS( SMESH::SMESH_IDSource_ptr meshPart, + const char* file, + CORBA::Boolean overwrite ) +{ +#ifdef WITH_CGNS + SMESH_TRY; + if ( _preMeshInfo ) + _preMeshInfo->FullLoadFromFile(); + + PrepareForWriting(file,overwrite); + + std::string meshName(""); + SALOMEDS::SObject_wrap so = _gen_i->ObjectToSObject( meshPart ); + if ( !so->_is_nil() ) + { + CORBA::String_var name = so->GetName(); + meshName = name.in(); + } + SMESH_TRY; + + SMESH::SMESH_Mesh_var mesh = meshPart->GetMesh(); + SMESH_Mesh_i* mesh_i = SMESH::DownCast( mesh ); + mesh_i->Load(); + auto myMesh = mesh_i->GetImpl().GetMeshDS(); + _impl->ExportStructuredCGNS(file, myMesh, meshName.c_str()); + + SMESH_CATCH( SMESH::throwCorbaException ); + + TPythonDump() << SMESH::SMESH_Mesh_var(_this()) << ".ExportStructuredCGNS( " + "r'" << file << "', " << overwrite << ", " << meshPart << ")"; + + SMESH_CATCH( SMESH::throwCorbaException ); + +#else + THROW_SALOME_CORBA_EXCEPTION("CGNS library is unavailable", SALOME::INTERNAL_ERROR); +#endif +} + //================================================================================ /*! * \brief Export a part of mesh to a GMF file diff --git a/src/SMESH_I/SMESH_Mesh_i.hxx b/src/SMESH_I/SMESH_Mesh_i.hxx index 131993a40..a7cb38fab 100644 --- a/src/SMESH_I/SMESH_Mesh_i.hxx +++ b/src/SMESH_I/SMESH_Mesh_i.hxx @@ -218,6 +218,11 @@ public: void ExportDAT( const char* file, const CORBA::Boolean renumber ); void ExportUNV( const char* file, const CORBA::Boolean renumber ); void ExportSTL( const char* file, bool isascii ); + + void ExportStructuredCGNS(SMESH::SMESH_IDSource_ptr meshPart, + const char* file, + CORBA::Boolean overwrite); + void ExportCGNS(SMESH::SMESH_IDSource_ptr meshPart, const char* file, CORBA::Boolean overwrite, @@ -268,7 +273,7 @@ public: CORBA::Boolean renumber); void ExportPartToSTL(SMESH::SMESH_IDSource_ptr meshPart, const char* file, - CORBA::Boolean isascii); + CORBA::Boolean isascii); CORBA::Double GetComputeProgress(); diff --git a/src/SMESH_SWIG/smeshBuilder.py b/src/SMESH_SWIG/smeshBuilder.py index ce1bfdc13..607314411 100644 --- a/src/SMESH_SWIG/smeshBuilder.py +++ b/src/SMESH_SWIG/smeshBuilder.py @@ -2588,6 +2588,26 @@ class Mesh(metaclass = MeshMeta): meshPart = self.mesh self.mesh.ExportCGNS(meshPart, f, overwrite, groupElemsByType) + def ExportStructuredCGNS(self, f, overwrite=1, meshPart=None ): + """ + Export the mesh in a file in Structured CGNS format + + Parameters: + f: is the file name + overwrite: boolean parameter for overwriting/not overwriting the file + meshPart: a part of mesh (:class:`sub-mesh, group or filter `) to export instead of the mesh + """ + + unRegister = genObjUnRegister() + if isinstance( meshPart, list ): + meshPart = self.GetIDSource( meshPart, SMESH.ALL ) + unRegister.set( meshPart ) + if isinstance( meshPart, Mesh ): + meshPart = meshPart.mesh + elif not meshPart: + meshPart = self.mesh + self.mesh.ExportStructuredCGNS(meshPart, f, overwrite) + def ExportGMF(self, f, meshPart=None): """ Export the mesh in a file in GMF format. diff --git a/src/StdMeshers/StdMeshers_Hexa_3D.cxx b/src/StdMeshers/StdMeshers_Hexa_3D.cxx index f7e43837f..8020d349f 100644 --- a/src/StdMeshers/StdMeshers_Hexa_3D.cxx +++ b/src/StdMeshers/StdMeshers_Hexa_3D.cxx @@ -897,7 +897,7 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, // compute internal node coordinates gp_XYZ coords; SMESH_Block::ShellPoint( params, pointsOnShapes, coords ); - column[ z ] = helper.AddNode( coords.X(), coords.Y(), coords.Z() ); + column[ z ] = helper.AddNode( coords.X(), coords.Y(), coords.Z() ); } // z loop if ( toRenumber ) @@ -926,7 +926,6 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, // 5) Create hexahedrons // --------------------- - for ( x = 0; x < xSize-1; ++x ) { for ( y = 0; y < ySize-1; ++y ) { vector< const SMDS_MeshNode* >& col00 = columns[ colIndex( x, y )]; @@ -936,12 +935,16 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, for ( z = 0; z < zSize-1; ++z ) { // bottom face normal of a hexa mush point outside the volume - if ( toRenumber ) + if ( toRenumber ) + { helper.AddVolume(col00[z], col01[z], col01[z+1], col00[z+1], - col10[z], col11[z], col11[z+1], col10[z+1]); + col10[z], col11[z], col11[z+1], col10[z+1]); + } else + { helper.AddVolume(col00[z], col01[z], col11[z], col10[z], - col00[z+1], col01[z+1], col11[z+1], col10[z+1]); + col00[z+1], col01[z+1], col11[z+1], col10[z+1]); + } } } } @@ -949,6 +952,11 @@ bool StdMeshers_Hexa_3D::Compute(SMESH_Mesh & aMesh, if ( toRenumber ) renumHelper.DoReplaceNodes(); + meshDS->SetStructuredGrid( aShape, zSize, ySize, xSize ); + for ( x = 0; x < xSize; ++x ) + for ( y = 0; y < ySize; ++y ) + for ( z = 0; z < zSize; ++z ) + meshDS->SetNodeOnStructuredGrid( aShape, columns[colIndex(x,y)][z], z, y, x ); if ( _blockRenumberHyp ) { diff --git a/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx b/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx index 1abd26595..7c27720d1 100644 --- a/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx +++ b/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx @@ -431,6 +431,16 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, SMESHDS_Mesh * meshDS = aMesh.GetMeshDS(); Handle(Geom_Surface) S = BRep_Tool::Surface(aFace); int i,j, geomFaceID = meshDS->ShapeToIndex(aFace); + + meshDS->SetStructuredGrid( aFace, nbhoriz, nbvertic ); + for (j = 0; j < nbvertic; j++) + for (i = 0; i < nbhoriz; i++) + { + UVPtStruct& uvPnt = quad->UVPt( i, j ); + auto P = std::make_shared( S->Value( uvPnt.u, uvPnt.v ).Coord() ); + meshDS->SetNodeOnStructuredGrid( aFace, P, i, j ); + } + for (i = 1; i < nbhoriz - 1; i++) for (j = 1; j < nbvertic - 1; j++) { diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 70a9e0e57..d52fdcf89 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -46,3 +46,22 @@ FOREACH(_test ${GOOD_TESTS}) COMMAND ${PYTHON_EXECUTABLE} -B ${CMAKE_SOURCE_DIR}/test/test_helper.py ${CMAKE_CURRENT_SOURCE_DIR}/${_test}) SET_TESTS_PROPERTIES(${testname} PROPERTIES ENVIRONMENT "${tests_env}" LABELS "tests") ENDFOREACH() + +INCLUDE_DIRECTORIES( + ${PROJECT_SOURCE_DIR}/src/SMESHUtils + ${PROJECT_SOURCE_DIR}/src/SMDS + ) + +FOREACH(_test ${CPP_TESTS}) + GET_FILENAME_COMPONENT(testname ${_test} NAME_WE) + SET(testname "TESTS_${testname}") + + add_executable(${_test} ${_test}.cxx) + target_link_libraries(${_test} SMESHUtils SMDS ) + + ADD_TEST(NAME ${testname} + COMMAND ${CMAKE_CURRENT_SOURCE_DIR}/${_test} ) + SET_TESTS_PROPERTIES(${testname} PROPERTIES ENVIRONMENT "${tests_env}" LABELS "tests") +ENDFOREACH() + +INSTALL(FILES ${CMAKE_CURRENT_BINARY_DIR}/${CPP_TESTS} PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ DESTINATION ${TEST_INSTALL_DIRECTORY}) diff --git a/test/CTestTestfileInstall.cmake b/test/CTestTestfileInstall.cmake index 4961f25d5..f4f1fed9f 100644 --- a/test/CTestTestfileInstall.cmake +++ b/test/CTestTestfileInstall.cmake @@ -34,3 +34,10 @@ FOREACH(tfile ${SESSION_FREE_TESTS}) ADD_TEST(${TEST_NAME} python ${tfile}) SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "${COMPONENT_NAME};${COMPONENT_NAME}_tests") endforeach() + +FOREACH(tfile ${CPP_TESTS}) + GET_FILENAME_COMPONENT(BASE_NAME ${tfile} NAME_WE) + SET(TEST_NAME SMESH_${BASE_NAME}) + ADD_TEST(${TEST_NAME} SMESH_RegularGridTest ) + SET_TESTS_PROPERTIES(${TEST_NAME} PROPERTIES LABELS "${COMPONENT_NAME};${COMPONENT_NAME}_tests") +endforeach() diff --git a/test/SMESH_RegularGridTest.cxx b/test/SMESH_RegularGridTest.cxx new file mode 100644 index 000000000..8d3fd930f --- /dev/null +++ b/test/SMESH_RegularGridTest.cxx @@ -0,0 +1,326 @@ +// Copyright (C) 2016-2024 CEA, EDF +// +// 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_RegularGridTest.cxx (unit test) +// Author : Cesar Conopoima (cce) + +// std +#include +#include + +// smesh +#include "SMESH_RegularGridTemplate.hxx" +#include "SMDS_MeshNode.hxx" + +bool testConstructor() +{ + auto regularGrid2D = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,5,4)); + if ( regularGrid2D->Size() != 5*4 ) throw std::runtime_error("2D Grid size not match in testConstructor()\n"); + auto regularGrid3D = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,5,4,20)); + if ( regularGrid3D->Size() != 5*4*20 ) throw std::runtime_error("3D Grid size not match in testConstructor()\n"); + return true; +} + +bool testSetGetNode() +{ + auto regularGrid = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,5,4)); + double x,y; + + for (size_t j = 0; j < 4; j++) + for (size_t i = 0; i < 5; i++) + { + x = i * 0.1; + y = j * 0.1; + auto node = std::make_shared( x, y, 0. ); + regularGrid->SetNode( node, i, j ); + } + for (size_t j = 0; j < 4; j++) + for (size_t i = 0; i < 5; i++) + { + x = i * 0.1; + y = j * 0.1; + auto node = regularGrid->GetNode( i, j, 0 ); + if ( node == nullptr ) throw std::runtime_error("error getting node in testSetGetNode()\n"); + if ( node->X() != x || + node->Y() != y || + node->Z() != 0. ) throw std::runtime_error("error for node not in the correct position in testSetGetNode()\n"); + } + + return true; +} + +bool testGetCommontInterface1D() +{ + int mnx = 5, mny = 4; + auto regularGrid0 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,mnx,mny)); + auto regularGrid1 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(1,mnx,mny)); + double x,y; + + /*side by side grid with 1to1 interface*/ + /* o-o-o-o-oxo-o-o-o-o */ + /* o-oG0-o-oxo-o-G1o-o */ + /* o-o-o-o-oxo-o-o-o-o */ + /* o-o-o-o-oxo-o-o-o-o */ + + for (size_t j = 0; j < mny; j++) + for (size_t i = 0; i < mnx; i++) + { + x = i * 0.1; + y = j * 0.1; + auto node0 = std::make_shared(x,y,0.); + regularGrid0->SetNode( node0, i, j ); + + x = i * 0.1 + 0.4; + y = j * 0.1; + auto node1 = std::make_shared(x,y,0.); + regularGrid1->SetNode( node1, i, j ); + } + + std::vector interface; /* interface at right */ + regularGrid0->GetEdgeInterfaces( regularGrid1.get(), interface ); + std::vector expectedRange = regularGrid0->getEdgeIndexLimits(SMESHUtils::SMESH_RegularGrid::RIGHT); + std::vector expectedDonorRange = regularGrid1->getEdgeIndexLimitsInverted(SMESHUtils::SMESH_RegularGrid::LEFT); + + if ( interface.size() != 1+mny*2+2 ) + throw std::runtime_error("1to1 edge interface not found in testGetCommontInterface1D()\n"); + + if ( interface[0] != int(SMESHUtils::SMESH_RegularGrid::RIGHT) ) + throw std::runtime_error("edge found in interface is incorrect testGetCommontInterface1D()\n"); + + for (size_t i = 0; i < 4; i++) + if ( interface[i+1] != expectedRange[i]) + throw std::runtime_error("interface range index is incorrect in testGetCommontInterface1D()\n"); + + for (size_t i = 0; i < 4; i++) + if ( interface[i+5] != expectedDonorRange[i]) + throw std::runtime_error("interface donor index is incorrect in testGetCommontInterface1D()\n"); + + return true; +} + +bool testGetCommontInterface2D() +{ + int mnx = 5, mny = 4, mnz = 2; + auto regularGrid0 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,mnx,mny,mnz)); + auto regularGrid1 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(1,mnx,mny,mnz)); + auto regularGrid2 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(2,mnx,mny,mnz)); + auto regularGrid3 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(3,mnx,mny,mnz)); + auto regularGrid4 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(4,mnx,mny,mnz)); + auto regularGrid5 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(5,mnx,mny,mnz)); + auto regularGrid6 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(6,mnx,mny,mnz)); + double x,y,z; + + for (size_t k = 0; k < mnz; k++) + for (size_t j = 0; j < mny; j++) + for (size_t i = 0; i < mnx; i++) + { + x = i * 0.1; + y = j * 0.1; + z = k * 0.1; + auto node0 = std::make_shared(x,y,z); + regularGrid0->SetNode( node0, i, j, k ); + + // Grid at botton + x = i * 0.1; + y = j * 0.1; + z = k * 0.1-0.1; + auto node1 = std::make_shared(x,y,z); + regularGrid1->SetNode( node1, i, j, k ); + + // Grid at right + x = i * 0.1+0.4; + y = j * 0.1; + z = k * 0.1; + auto node2 = std::make_shared(x,y,z); + regularGrid2->SetNode( node2, i, j, k ); + + // Grid at back + x = i * 0.1; + y = j * 0.1+0.3; + z = k * 0.1; + auto node3 = std::make_shared(x,y,z); + regularGrid3->SetNode( node3, i, j, k ); + + // Grid at left + x = i * 0.1-0.4; + y = j * 0.1; + z = k * 0.1; + auto node4 = std::make_shared(x,y,z); + regularGrid4->SetNode( node4, i, j, k ); + + // Grid at front + x = i * 0.1; + y = j * 0.1-0.3; + z = k * 0.1; + auto node5 = std::make_shared(x,y,z); + regularGrid5->SetNode( node5, i, j, k ); + + // Grid on top + x = i * 0.1; + y = j * 0.1; + z = k * 0.1+0.1; + auto node6 = std::make_shared(x,y,z); + regularGrid6->SetNode( node6, i, j, k ); + } + + std::vector> interface(6); + std::vector> expectedRange(6); + std::vector> expectedDonorRange(6); + + regularGrid0->GetFaceInterfaces( regularGrid1.get(), interface[0] ); + regularGrid0->GetFaceInterfaces( regularGrid2.get(), interface[1] ); + regularGrid0->GetFaceInterfaces( regularGrid3.get(), interface[2] ); + regularGrid0->GetFaceInterfaces( regularGrid4.get(), interface[3] ); + regularGrid0->GetFaceInterfaces( regularGrid5.get(), interface[4] ); + regularGrid0->GetFaceInterfaces( regularGrid6.get(), interface[5] ); + + expectedRange[0] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_BOTTOM); + expectedRange[1] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_RIGHT); + expectedRange[2] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_BACK); + expectedRange[3] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_LEFT); + expectedRange[4] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_FRONT); + expectedRange[5] = regularGrid0->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_TOP); + + expectedDonorRange[0] = regularGrid1->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_TOP); + expectedDonorRange[1] = regularGrid2->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_LEFT); + expectedDonorRange[2] = regularGrid3->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_FRONT ); + expectedDonorRange[3] = regularGrid4->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_RIGHT); + expectedDonorRange[4] = regularGrid5->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_BACK); + expectedDonorRange[5] = regularGrid6->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_BOTTOM); + + int interfaceSize = 1+6*2+3; + for ( auto itf : interface ) + { + if ( itf.size() != interfaceSize ) + throw std::runtime_error("1to1 face interface not found in testGetCommontInterface2D()\n"); + } + + int count = 0; + for ( auto itf : interface ) + { + if ( itf[0] != count /*face id of the interface*/) + throw std::runtime_error("face found in interface is incorrect testGetCommontInterface2D()\n"); + + for (size_t i = 0; i < 6; i++) + if ( itf[i+1] != expectedRange[count][i]) + throw std::runtime_error("interface range index is incorrect in testGetCommontInterface2D()\n"); + + for (size_t i = 0; i < 6; i++) + if ( itf[i+7] != expectedDonorRange[count][i]) + throw std::runtime_error("interface donor index is incorrect in testGetCommontInterface2D()\n"); + count++; + } + + return true; +} + +bool testGetPartialInterface2D() +{ + int mnx = 5, mny = 4, mnz = 3; + int mnxf = 2, mnyf = 2, mnzf = 2; + int mnxl = 2, mnyl = 3, mnzl = 2; + auto regularGrid0 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(0,mnx,mny,mnz)); + auto regularGrid1 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(1,mnxf,mnyf,mnzf)); + auto regularGrid2 = std::unique_ptr( new SMESHUtils::SMESH_RegularGrid(2,mnxl,mnyl,mnzl)); + double x,y,z; + + for (size_t k = 0; k < mnz; k++) + for (size_t j = 0; j < mny; j++) + for (size_t i = 0; i < mnx; i++) + { + x = i * 0.1; + y = j * 0.1; + z = k * 0.1; + auto node0 = std::make_shared(x,y,z); + regularGrid0->SetNode( node0, i, j, k ); + } + + for (size_t k = 0; k < mnzf; k++) + for (size_t j = 0; j < mnyf; j++) + for (size_t i = 0; i < mnxf; i++) + { + // Grid at front + x = i * 0.1+0.3; + y = -j * 0.1; + z = k * 0.1+0.1; + auto node1 = std::make_shared(x,y,z); + regularGrid1->SetNode( node1, i, j, k ); + } + + for (size_t k = 0; k < mnzl; k++) + for (size_t j = 0; j < mnyl; j++) + for (size_t i = 0; i < mnxl; i++) + { + // Grid at left + x = -i * 0.1; + y = j * 0.1+0.1; + z = k * 0.1; + auto node2 = std::make_shared(x,y,z); + regularGrid2->SetNode( node2, i, j, k ); + } + + std::vector> interface(2); + std::vector faceId = {SMESHUtils::SMESH_RegularGrid::B_FRONT,SMESHUtils::SMESH_RegularGrid::B_LEFT}; /*faces where interfaces are*/ + std::vector> expectedRange(2); + std::vector> expectedDonorRange(2); + + regularGrid0->GetFaceInterfaces( regularGrid1.get(), interface[0] ); + regularGrid0->GetFaceInterfaces( regularGrid2.get(), interface[1] ); + + expectedRange[0] = std::vector({4,1,2,5,1,3}); + expectedRange[1] = std::vector({1,4,1,1,2,2}); + + expectedDonorRange[0] = regularGrid1->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_FRONT); + expectedDonorRange[1] = regularGrid2->getFaceIndexLimits(SMESHUtils::SMESH_RegularGrid::B_LEFT); + + int interfaceSize = 1+6*2+3; + for ( auto itf : interface ) + { + if ( itf.size() != interfaceSize ) + throw std::runtime_error("1to1 face interface not found in testGetPartialInterface2D()\n"); + } + + int count = 0; + for ( auto itf : interface ) + { + if ( itf[0] != faceId[count] /*face id of the interface*/) + throw std::runtime_error("face found in interface is incorrect testGetPartialInterface2D()\n"); + + for (size_t i = 0; i < 6; i++) + if ( itf[i+1] != expectedRange[count][i]) + throw std::runtime_error("interface range index is incorrect in testGetPartialInterface2D()\n"); + + for (size_t i = 0; i < 6; i++) + if ( itf[i+7] != expectedDonorRange[count][i]) + throw std::runtime_error("interface donor index is incorrect in testGetPartialInterface2D()\n"); + count++; + } + + return true; +} + +int main() +{ + if ( !testConstructor() || !testSetGetNode() || + !testGetCommontInterface1D() || !testGetCommontInterface2D() || + !testGetPartialInterface2D() ) + return 1; + else + return 0; + +} \ No newline at end of file diff --git a/test/tests.set b/test/tests.set index 9ed691711..f1a155efc 100644 --- a/test/tests.set +++ b/test/tests.set @@ -122,6 +122,11 @@ SET(GOOD_TESTS SMESH_test5.py ) + +SET(CPP_TESTS + SMESH_RegularGridTest +) + # The following tests can be executed without driver, just by python. # ---------------------------------------------------------------------------- -- 2.30.2