diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 367d0bb2662..5fa42d5fc4b 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -20,6 +20,8 @@ #include "common/DataLayouts.hpp" #include "common/TimingMacros.hpp" +#include "mesh/WellElementRegion.hpp" +#include "mesh/WellElementSubRegion.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" #include "SurfaceElementRegion.hpp" #include "constitutive/ConstitutiveManager.hpp" @@ -228,6 +230,33 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); + meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion const & region ){ + region.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion const & subRegion ){ + TableData dataPerforation; + TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", region.getName()), + { "Perforation no.", "Well element no.", "Coordinates", + "Cell region", "Cell sub-region", "Cell ID" } ); + for( globalIndex iperf = 0; iperf < subRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) + { + if( subRegion.getReservoirElementID() != -1 ) + { + arrayView2d< const real64 > const perfLocation = subRegion.getPerforationData()->getLocation(); + dataPerforation.addRow( iperf, + subRegion.getPerforationData()->getWellElements()[iperf], + perfLocation[iperf], + subRegion.getRegionName(), + subRegion.getSubRegionName(), subRegion.getReservoirElementID()); + + } + } + if( dataPerforation.getCellsData().size() > 0 ) + { + TableTextFormatter const formatter( layoutPerforation ); + GEOS_LOG( formatter.toString( dataPerforation )); + } + } ); + } ); + // communicate to rebuild global node info since we modified global ordering nodeManager.setMaxGlobalIndex(); } diff --git a/src/coreComponents/mesh/WellElementSubRegion.cpp b/src/coreComponents/mesh/WellElementSubRegion.cpp index b58343bdd12..66e9c4babaa 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.cpp +++ b/src/coreComponents/mesh/WellElementSubRegion.cpp @@ -33,7 +33,10 @@ WellElementSubRegion::WellElementSubRegion( string const & name, Group * const p m_topWellElementIndex( -1 ), m_perforationData( groupKeyStruct::perforationDataString(), this ), m_topRank( -1 ), - m_searchDepth( 10 ) + m_searchDepth( 10 ), + m_reservoirElementID(-1), + m_regionName( "" ), + m_subRegionName( "" ) { m_elementType = ElementType::Line; @@ -377,107 +380,6 @@ void initializeLocalSearch( MeshLevel const & mesh, // note that this reservoir element does not necessarily contains "location" eiInit = ret.second; } - -/** - * @brief Search for the reservoir element that contains the well element. - To do that, loop over the reservoir elements that are in the neighborhood of (erInit,esrInit,eiInit) - * @param[in] meshLevel the mesh object (single level only) - * @param[in] location the location of that we are trying to match with a reservoir element - * @param[in] targetRegionIndex the target region index for the reservoir element - * @param[inout] esrMatched the subregion index of the reservoir element that contains "location", if any - * @param[inout] eiMatched the element index of the reservoir element that contains "location", if any - * @param[inout] giMatched the element global index of the reservoir element that contains "location", if any - */ -bool searchLocalElements( MeshLevel const & mesh, - real64 const (&location)[3], - localIndex const & searchDepth, - localIndex const & targetRegionIndex, - localIndex & esrMatched, - localIndex & eiMatched, - globalIndex & giMatched, - real64 const geomTol ) -{ - ElementRegionBase const & region = mesh.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); - - bool resElemFound = false; - for( localIndex esr = 0; esr < region.numSubRegions(); ++esr ) - { - ElementSubRegionBase const & subRegionBase = region.getSubRegion( esr ); - region.applyLambdaToContainer< CellElementSubRegion, SurfaceElementSubRegion >( subRegionBase, [&]( auto const & subRegion ) - { - GEOS_LOG_RANK_0( GEOS_FMT( " searching well connections with region/subregion: {}/{}", region.getName(), subRegion.getName() ) ); - - // first, we search for the reservoir element that is the *closest* from the center of well element - // note that this reservoir element does not necessarily contain the center of the well element - // this "init" reservoir element will be used later to find the reservoir element that - // contains the well element - localIndex eiInit = -1; - - initializeLocalSearch( mesh, location, targetRegionIndex, esr, eiInit ); - - if( eiInit < 0 ) // nothing found, skip the rest - return; - - // loop over the reservoir elements that are in the neighborhood of (esrInit,eiInit) - // search locally, starting from the location of the previous perforation - // the assumption here is that perforations have been entered in order of depth - - SortedArray< localIndex > nodes; - SortedArray< globalIndex > elements; - - // here is how the search is done: - // 1 - We check if "location" is within the "init" reservoir element defined by (erInit,esrMatched,eiMatched) - // 2 - If yes, stop - // - If not, a) collect the nodes of the reservoir element defined by (erInit,esrMatched,eiMatched) - // b) use these nodes to grab the neighbors of (erInit,esrMatched,eiMatched) - // c) check if "location" is within the neighbors. If not, grab the neighbors of the neighbors, and so - // on... - - // collect the nodes of the current element - // they will be used to access the neighbors and check if they contain the perforation - collectElementNodes( subRegion, eiInit, nodes ); - - // if no match is found, enlarge the neighborhood m_searchDepth'th times - for( localIndex d = 0; d < searchDepth; ++d ) - { - localIndex nNodes = nodes.size(); - - // search the reservoir elements that can be accessed from the set "nodes" - // stop if a reservoir element containing the perforation is found - // if not, enlarge the set "nodes" - - resElemFound = - visitNeighborElements< TYPEOFREF( subRegion ) >( mesh, - location, - nodes, - elements, - targetRegionIndex, - esr, - eiMatched, - giMatched, - geomTol ); - - if( resElemFound || nNodes == nodes.size()) - { - if( resElemFound ) - { - esrMatched = esr; - GEOS_LOG( GEOS_FMT( " found {}/{}/{}", region.getName(), subRegion.getName(), giMatched ) ); - } - return; - } - } - } ); - - if( resElemFound ) - { - break; - } - } - - return resElemFound; -} - } void WellElementSubRegion::generate( MeshLevel & mesh, @@ -583,12 +485,106 @@ void WellElementSubRegion::generate( MeshLevel & mesh, } +bool WellElementSubRegion::searchLocalElements( MeshLevel const & mesh, + real64 const (&location)[3], + localIndex const & searchDepth, + localIndex const & targetRegionIndex, + localIndex & esrMatched, + localIndex & eiMatched, + globalIndex & giMatched, + real64 const geomTol ) +{ + ElementRegionBase const & region = mesh.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); + + bool resElemFound = false; + for( localIndex esr = 0; esr < region.numSubRegions(); ++esr ) + { + ElementSubRegionBase const & subRegionBase = region.getSubRegion( esr ); + region.applyLambdaToContainer< CellElementSubRegion, SurfaceElementSubRegion >( subRegionBase, [&]( auto const & subRegion ) + { + GEOS_LOG_RANK_0( GEOS_FMT( " searching well connections with region/subregion: {}/{}", region.getName(), subRegion.getName() ) ); + + // first, we search for the reservoir element that is the *closest* from the center of well element + // note that this reservoir element does not necessarily contain the center of the well element + // this "init" reservoir element will be used later to find the reservoir element that + // contains the well element + localIndex eiInit = -1; + + initializeLocalSearch( mesh, location, targetRegionIndex, esr, eiInit ); + + if( eiInit < 0 ) // nothing found, skip the rest + return; + + // loop over the reservoir elements that are in the neighborhood of (esrInit,eiInit) + // search locally, starting from the location of the previous perforation + // the assumption here is that perforations have been entered in order of depth + + SortedArray< localIndex > nodes; + SortedArray< globalIndex > elements; + + // here is how the search is done: + // 1 - We check if "location" is within the "init" reservoir element defined by (erInit,esrMatched,eiMatched) + // 2 - If yes, stop + // - If not, a) collect the nodes of the reservoir element defined by (erInit,esrMatched,eiMatched) + // b) use these nodes to grab the neighbors of (erInit,esrMatched,eiMatched) + // c) check if "location" is within the neighbors. If not, grab the neighbors of the neighbors, and so + // on... + + // collect the nodes of the current element + // they will be used to access the neighbors and check if they contain the perforation + collectElementNodes( subRegion, eiInit, nodes ); + + // if no match is found, enlarge the neighborhood m_searchDepth'th times + for( localIndex d = 0; d < searchDepth; ++d ) + { + localIndex nNodes = nodes.size(); + + // search the reservoir elements that can be accessed from the set "nodes" + // stop if a reservoir element containing the perforation is found + // if not, enlarge the set "nodes" + + resElemFound = + visitNeighborElements< TYPEOFREF( subRegion ) >( mesh, + location, + nodes, + elements, + targetRegionIndex, + esr, + eiMatched, + giMatched, + geomTol ); + + if( resElemFound || nNodes == nodes.size()) + { + if( resElemFound ) + { + esrMatched = esr; + m_reservoirElementID = giMatched; + m_regionName = region.getName(); + m_subRegionName = subRegion.getName(); + GEOS_LOG( GEOS_FMT( " found {}/{}/{}", region.getName(), subRegion.getName(), giMatched ) ); + } + return; + } + } + } ); + + if( resElemFound ) + { + break; + } + } + + return resElemFound; +} + + void WellElementSubRegion::assignUnownedElementsInReservoir( MeshLevel & mesh, LineBlockABC const & lineBlock, SortedArray< globalIndex > const & unownedElems, SortedArray< globalIndex > & localElems, arrayView1d< integer > & elemStatusGlobal, - real64 const geomTol ) const + real64 const geomTol ) { ElementRegionManager const & elemManager = mesh.getElemManager(); // get the well and reservoir element coordinates diff --git a/src/coreComponents/mesh/WellElementSubRegion.hpp b/src/coreComponents/mesh/WellElementSubRegion.hpp index 26089a42838..addae45e5f2 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.hpp +++ b/src/coreComponents/mesh/WellElementSubRegion.hpp @@ -383,6 +383,25 @@ class WellElementSubRegion : public ElementSubRegionBase * @return list of indicies */ array1d< globalIndex > const & getGlobalElementIndex() const { return m_globalElementIndex; } + + /** + * @return The reservoir element that contains the perforation + */ + globalIndex getReservoirElementID() const + { return m_reservoirElementID; } + + /** + * @return The region name containing the reservoir element + */ + string getSubRegionName() const + { return m_subRegionName; } + + /** + * @return The sub region name containing the reservoir element + */ + string getRegionName() const + { return m_regionName; } + private: /** @@ -403,7 +422,7 @@ class WellElementSubRegion : public ElementSubRegionBase SortedArray< globalIndex > const & unownedElems, SortedArray< globalIndex > & localElems, arrayView1d< integer > & elemStatusGlobal, - real64 geomTol ) const; + real64 geomTol ); /** * @brief Check that all the well elements have been assigned to a single rank. @@ -461,6 +480,26 @@ class WellElementSubRegion : public ElementSubRegionBase */ void updateNodeManagerNodeToElementMap( MeshLevel & mesh ); + /** + * @brief Search for the reservoir element that contains the well element. + To do that, loop over the reservoir elements that are in the neighborhood of (erInit,esrInit,eiInit) + * @param[in] meshLevel the mesh object (single level only) + * @param[in] location the location of that we are trying to match with a reservoir element + * @param[in] targetRegionIndex the target region index for the reservoir element + * @param[inout] esrMatched the subregion index of the reservoir element that contains "location", if any + * @param[inout] eiMatched the element index of the reservoir element that contains "location", if any + * @param[inout] giMatched the element global index of the reservoir element that contains "location", if any + */ + bool searchLocalElements( MeshLevel const & mesh, + real64 const (&location)[3], + localIndex const & searchDepth, + localIndex const & targetRegionIndex, + localIndex & esrMatched, + localIndex & eiMatched, + globalIndex & giMatched, + real64 const geomTol ); + + /** * @brief Pack element-to-node and element-to-face maps * @tparam the flag for the bufferOps::Pack function @@ -522,6 +561,15 @@ class WellElementSubRegion : public ElementSubRegionBase /// Number of segment per rank array1d< localIndex > m_elementPerRank; + + /// The reservoir element that contains the perforation + globalIndex m_reservoirElementID; + /// the target region name for the reservoir element + string m_regionName; + /// the target sub region name for the reservoir element + string m_subRegionName; + + }; } /* namespace geos */ diff --git a/src/coreComponents/mesh/generators/WellGeneratorBase.cpp b/src/coreComponents/mesh/generators/WellGeneratorBase.cpp index f55ba78049d..2f8c37a7b12 100644 --- a/src/coreComponents/mesh/generators/WellGeneratorBase.cpp +++ b/src/coreComponents/mesh/generators/WellGeneratorBase.cpp @@ -140,7 +140,6 @@ void WellGeneratorBase::generateWellGeometry( ) if( isLogLevelActive< logInfo::GenerateWell >( this->getLogLevel() ) && MpiWrapper::commRank() == 0 ) { logInternalWell(); - logPerforationTable(); } } @@ -557,18 +556,4 @@ void WellGeneratorBase::logInternalWell() const GEOS_LOG_RANK_0( wellFormatter.toString( wellData )); } -void WellGeneratorBase::logPerforationTable() const -{ - TableData dataPerforation; - for( globalIndex iperf = 0; iperf < m_numPerforations; ++iperf ) - { - dataPerforation.addRow( iperf, m_perfCoords[iperf], m_perfElemId[iperf] ); - } - - TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", getName()), - { "Perforation no.", "Coordinates", "Well element no." } ); - TableTextFormatter const logPerforation( layoutPerforation ); - GEOS_LOG_RANK_0( logPerforation.toString( dataPerforation )); -} - } diff --git a/src/coreComponents/mesh/generators/WellGeneratorBase.hpp b/src/coreComponents/mesh/generators/WellGeneratorBase.hpp index e9a2dd81046..0367b2f0109 100644 --- a/src/coreComponents/mesh/generators/WellGeneratorBase.hpp +++ b/src/coreComponents/mesh/generators/WellGeneratorBase.hpp @@ -306,7 +306,6 @@ class WellGeneratorBase : public MeshComponentBase /// @cond DO_NOT_DOCUMENT void logInternalWell() const; - void logPerforationTable() const; /// @endcond /// Global number of perforations