From 87a4eecfdf9f3b626d8836221d8334358f172f5c Mon Sep 17 00:00:00 2001 From: arng40 Date: Fri, 6 Feb 2026 11:41:18 +0100 Subject: [PATCH 1/8] first attempt new version of log table perforation --- .../mesh/ElementRegionManager.cpp | 26 +++ .../mesh/WellElementSubRegion.cpp | 197 +++++++++--------- .../mesh/WellElementSubRegion.hpp | 51 ++++- .../mesh/generators/WellGeneratorBase.cpp | 15 -- 4 files changed, 171 insertions(+), 118 deletions(-) diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 367d0bb2662..055cae2fc76 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -20,6 +20,9 @@ #include "common/DataLayouts.hpp" #include "common/TimingMacros.hpp" +#include "common/logger/Logger.hpp" +#include "mesh/WellElementRegion.hpp" +#include "mesh/WellElementSubRegion.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" #include "SurfaceElementRegion.hpp" #include "constitutive/ConstitutiveManager.hpp" @@ -228,6 +231,29 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); + + meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion & region ){ + TableData dataPerforation; + TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", region.getName()), + { "Perforation no.", "Coordinates", "Well element no.", "Cell region", "Cell sub-region", "Cell ID" } ); + region.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion & subRegion ){ + for( globalIndex iperf = 0; iperf < subRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) + { + if( subRegion.getGlobalIndex() != -1 ) + { + arrayView2d< real64 > const perfLocation = subRegion.getPerforationData()->getLocation(); + dataPerforation.addRow( iperf, + perfLocation[iperf], + subRegion.getPerforationData()->getWellElements()[iperf], + subRegion.getRegionName(), + subRegion.getSubRegionName(), subRegion.getGlobalIndex()); + } + } + TableTextFormatter formatter( layoutPerforation ); + GEOS_LOG_RANK_0( 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..fbe888d168b 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.cpp +++ b/src/coreComponents/mesh/WellElementSubRegion.cpp @@ -377,107 +377,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 +482,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_perfoGlobalIndex = 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..00f67e8ae9a 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.hpp +++ b/src/coreComponents/mesh/WellElementSubRegion.hpp @@ -383,6 +383,29 @@ class WellElementSubRegion : public ElementSubRegionBase * @return list of indicies */ array1d< globalIndex > const & getGlobalElementIndex() const { return m_globalElementIndex; } + + /** + * @return The global index perfo ?? + */ + globalIndex getGlobalIndex() const + { + return m_perfoGlobalIndex; + } + /** + * @return todo + */ + string getSubRegionName() const + { + return m_subRegionName; + } + /** + * @return todo + */ + string getRegionName() const + { + return m_regionName; + } + private: /** @@ -403,7 +426,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 +484,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 +565,12 @@ class WellElementSubRegion : public ElementSubRegionBase /// Number of segment per rank array1d< localIndex > m_elementPerRank; + + globalIndex m_perfoGlobalIndex = -1; + string m_regionName; + 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 )); -} - } From 56380622b1d6fc20b8215c81defad6d41d669257 Mon Sep 17 00:00:00 2001 From: arng40 Date: Fri, 6 Feb 2026 16:11:20 +0100 Subject: [PATCH 2/8] doxygen + add const --- .../mesh/ElementRegionManager.cpp | 29 ++++---- .../mesh/WellElementSubRegion.cpp | 7 +- .../mesh/WellElementSubRegion.hpp | 69 +++++++++---------- 3 files changed, 55 insertions(+), 50 deletions(-) diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 055cae2fc76..5fa42d5fc4b 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -20,7 +20,6 @@ #include "common/DataLayouts.hpp" #include "common/TimingMacros.hpp" -#include "common/logger/Logger.hpp" #include "mesh/WellElementRegion.hpp" #include "mesh/WellElementSubRegion.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" @@ -231,26 +230,30 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); - - meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion & region ){ - TableData dataPerforation; - TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", region.getName()), - { "Perforation no.", "Coordinates", "Well element no.", "Cell region", "Cell sub-region", "Cell ID" } ); - region.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion & subRegion ){ + 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.getGlobalIndex() != -1 ) + if( subRegion.getReservoirElementID() != -1 ) { - arrayView2d< real64 > const perfLocation = subRegion.getPerforationData()->getLocation(); + arrayView2d< const real64 > const perfLocation = subRegion.getPerforationData()->getLocation(); dataPerforation.addRow( iperf, - perfLocation[iperf], subRegion.getPerforationData()->getWellElements()[iperf], + perfLocation[iperf], subRegion.getRegionName(), - subRegion.getSubRegionName(), subRegion.getGlobalIndex()); + subRegion.getSubRegionName(), subRegion.getReservoirElementID()); + } } - TableTextFormatter formatter( layoutPerforation ); - GEOS_LOG_RANK_0( formatter.toString( dataPerforation )); + if( dataPerforation.getCellsData().size() > 0 ) + { + TableTextFormatter const formatter( layoutPerforation ); + GEOS_LOG( formatter.toString( dataPerforation )); + } } ); } ); diff --git a/src/coreComponents/mesh/WellElementSubRegion.cpp b/src/coreComponents/mesh/WellElementSubRegion.cpp index fbe888d168b..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; @@ -556,7 +559,7 @@ bool WellElementSubRegion::searchLocalElements( MeshLevel const & mesh, if( resElemFound ) { esrMatched = esr; - m_perfoGlobalIndex = giMatched; + m_reservoirElementID = giMatched; m_regionName = region.getName(); m_subRegionName = subRegion.getName(); GEOS_LOG( GEOS_FMT( " found {}/{}/{}", region.getName(), subRegion.getName(), giMatched ) ); diff --git a/src/coreComponents/mesh/WellElementSubRegion.hpp b/src/coreComponents/mesh/WellElementSubRegion.hpp index 00f67e8ae9a..addae45e5f2 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.hpp +++ b/src/coreComponents/mesh/WellElementSubRegion.hpp @@ -385,26 +385,22 @@ class WellElementSubRegion : public ElementSubRegionBase array1d< globalIndex > const & getGlobalElementIndex() const { return m_globalElementIndex; } /** - * @return The global index perfo ?? - */ - globalIndex getGlobalIndex() const - { - return m_perfoGlobalIndex; - } + * @return The reservoir element that contains the perforation + */ + globalIndex getReservoirElementID() const + { return m_reservoirElementID; } + /** - * @return todo - */ + * @return The region name containing the reservoir element + */ string getSubRegionName() const - { - return m_subRegionName; - } + { return m_subRegionName; } + /** - * @return todo - */ + * @return The sub region name containing the reservoir element + */ string getRegionName() const - { - return m_regionName; - } + { return m_regionName; } private: @@ -485,24 +481,24 @@ class WellElementSubRegion : public ElementSubRegionBase void updateNodeManagerNodeToElementMap( MeshLevel & mesh ); /** - * @brief Search for the reservoir element that contains the well element. + * @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 ); - + * @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 @@ -566,11 +562,14 @@ bool searchLocalElements( MeshLevel const & mesh, /// Number of segment per rank array1d< localIndex > m_elementPerRank; - globalIndex m_perfoGlobalIndex = -1; + /// 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 */ From 3b7ade4ac515c786ec8021c29eeafb74ccf213dc Mon Sep 17 00:00:00 2001 From: arng40 Date: Fri, 6 Feb 2026 16:18:23 +0100 Subject: [PATCH 3/8] remove logPerfo header --- src/coreComponents/mesh/generators/WellGeneratorBase.hpp | 1 - 1 file changed, 1 deletion(-) 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 From dbee251be6abdb4e26819b96ee0ba4039c65f8a0 Mon Sep 17 00:00:00 2001 From: arng40 Date: Tue, 10 Feb 2026 11:46:27 +0100 Subject: [PATCH 4/8] revert modif on wellElementSubRegion --- .../mesh/WellElementSubRegion.cpp | 202 +++++++++--------- .../mesh/WellElementSubRegion.hpp | 50 +---- 2 files changed, 104 insertions(+), 148 deletions(-) diff --git a/src/coreComponents/mesh/WellElementSubRegion.cpp b/src/coreComponents/mesh/WellElementSubRegion.cpp index 66e9c4babaa..b58343bdd12 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.cpp +++ b/src/coreComponents/mesh/WellElementSubRegion.cpp @@ -33,10 +33,7 @@ WellElementSubRegion::WellElementSubRegion( string const & name, Group * const p m_topWellElementIndex( -1 ), m_perforationData( groupKeyStruct::perforationDataString(), this ), m_topRank( -1 ), - m_searchDepth( 10 ), - m_reservoirElementID(-1), - m_regionName( "" ), - m_subRegionName( "" ) + m_searchDepth( 10 ) { m_elementType = ElementType::Line; @@ -380,6 +377,107 @@ 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, @@ -485,106 +583,12 @@ 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 ) + real64 const geomTol ) const { 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 addae45e5f2..26089a42838 100644 --- a/src/coreComponents/mesh/WellElementSubRegion.hpp +++ b/src/coreComponents/mesh/WellElementSubRegion.hpp @@ -383,25 +383,6 @@ 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: /** @@ -422,7 +403,7 @@ class WellElementSubRegion : public ElementSubRegionBase SortedArray< globalIndex > const & unownedElems, SortedArray< globalIndex > & localElems, arrayView1d< integer > & elemStatusGlobal, - real64 geomTol ); + real64 geomTol ) const; /** * @brief Check that all the well elements have been assigned to a single rank. @@ -480,26 +461,6 @@ 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 @@ -561,15 +522,6 @@ 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 */ From b6f5ca6b59529990ea304fd0e4c96f639f6204c3 Mon Sep 17 00:00:00 2001 From: arng40 Date: Tue, 10 Feb 2026 11:47:03 +0100 Subject: [PATCH 5/8] new method to get cell id & region / subregion --- .../mesh/ElementRegionManager.cpp | 26 ++++++++++++------- 1 file changed, 17 insertions(+), 9 deletions(-) diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 5fa42d5fc4b..932812fc8f3 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -20,6 +20,7 @@ #include "common/DataLayouts.hpp" #include "common/TimingMacros.hpp" +#include "mesh/ElementSubRegionBase.hpp" #include "mesh/WellElementRegion.hpp" #include "mesh/WellElementSubRegion.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" @@ -230,22 +231,29 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); - meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion const & region ){ - region.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion const & subRegion ){ + meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion const & wellRegion ){ + wellRegion.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion const & wellSubRegion ){ TableData dataPerforation; - TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", region.getName()), + TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", wellRegion.getName()), { "Perforation no.", "Well element no.", "Coordinates", "Cell region", "Cell sub-region", "Cell ID" } ); - for( globalIndex iperf = 0; iperf < subRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) + for( globalIndex iperf = 0; iperf < wellSubRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) { - if( subRegion.getReservoirElementID() != -1 ) + globalIndex cellId = wellSubRegion.getPerforationData()->getReservoirElementGlobalIndex()[iperf]; + if( cellId != -1 ) { - arrayView2d< const real64 > const perfLocation = subRegion.getPerforationData()->getLocation(); + localIndex targetRegionIndex = wellSubRegion.getPerforationData()->getMeshElements().m_toElementRegion[iperf]; + localIndex targetSubRegionIndex = wellSubRegion.getPerforationData()->getMeshElements().m_toElementSubRegion[iperf]; + ElementRegionBase const & region = meshLevel.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); + ElementSubRegionBase const & subRegion = region.getSubRegion< ElementSubRegionBase >( targetSubRegionIndex ); + + arrayView2d< const real64 > const perfLocation = wellSubRegion.getPerforationData()->getLocation(); dataPerforation.addRow( iperf, - subRegion.getPerforationData()->getWellElements()[iperf], + wellSubRegion.getPerforationData()->getWellElements()[iperf], perfLocation[iperf], - subRegion.getRegionName(), - subRegion.getSubRegionName(), subRegion.getReservoirElementID()); + region.getName(), + subRegion.getName(), + cellId ); } } From f10caa55994fa2ab8315f1242f35ffd1bfc65d9a Mon Sep 17 00:00:00 2001 From: arng40 Date: Tue, 10 Feb 2026 17:38:54 +0100 Subject: [PATCH 6/8] wip MPI for table perfo --- .../mesh/ElementRegionManager.cpp | 72 +++++++++++++------ 1 file changed, 50 insertions(+), 22 deletions(-) diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 932812fc8f3..3011115ed4e 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -19,7 +19,12 @@ #include "ElementRegionManager.hpp" #include "common/DataLayouts.hpp" +#include "common/DataTypes.hpp" +#include "common/MpiWrapper.hpp" +#include "common/StdContainerWrappers.hpp" #include "common/TimingMacros.hpp" +#include "common/logger/Logger.hpp" +#include "dataRepository/BufferOps.hpp" #include "mesh/ElementSubRegionBase.hpp" #include "mesh/WellElementRegion.hpp" #include "mesh/WellElementSubRegion.hpp" @@ -32,6 +37,8 @@ #include "schema/schemaUtilities.hpp" #include "mesh/generators/LineBlockABC.hpp" #include "mesh/CellElementRegionSelector.hpp" +#include "dataRepository/BufferOps.hpp" +#include "dataRepository/BufferOps_inline.hpp" namespace geos @@ -231,38 +238,59 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); - meshLevel.getElemManager().forElementRegions< WellElementRegion >( [&]( WellElementRegion const & wellRegion ){ - wellRegion.forElementSubRegions< WellElementSubRegion >( [&]( WellElementSubRegion const & wellSubRegion ){ - TableData dataPerforation; - TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", wellRegion.getName()), - { "Perforation no.", "Well element no.", "Coordinates", - "Cell region", "Cell sub-region", "Cell ID" } ); - for( globalIndex iperf = 0; iperf < wellSubRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) + forElementRegions< WellElementRegion >( [&]( WellElementRegion const & wellRegion ){ + + WellElementSubRegion const & + wellSubRegion = wellRegion.getGroup( ElementRegionBase::viewKeyStruct::elementSubRegions() ) + .getGroup< WellElementSubRegion >( wellRegion.getSubRegionName() ); + TableData dataPerforation; + TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", wellRegion.getWellGeneratorName()), + { "Perforation no.", "Well element no.", "Coordinates", + "Cell region", "Cell sub-region", "Cell ID" } ); + for( globalIndex iperf = 0; iperf < wellSubRegion.getPerforationData()->getNumPerforationsGlobal(); ++iperf ) + { + globalIndex const cellId = wellSubRegion.getPerforationData()->getReservoirElementGlobalIndex()[iperf]; + if( cellId != -1 ) { - globalIndex cellId = wellSubRegion.getPerforationData()->getReservoirElementGlobalIndex()[iperf]; - if( cellId != -1 ) - { - localIndex targetRegionIndex = wellSubRegion.getPerforationData()->getMeshElements().m_toElementRegion[iperf]; - localIndex targetSubRegionIndex = wellSubRegion.getPerforationData()->getMeshElements().m_toElementSubRegion[iperf]; - ElementRegionBase const & region = meshLevel.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); - ElementSubRegionBase const & subRegion = region.getSubRegion< ElementSubRegionBase >( targetSubRegionIndex ); + auto const & meshElems = wellSubRegion.getPerforationData()->getMeshElements(); + localIndex const targetRegionIndex = meshElems.m_toElementRegion[iperf]; + localIndex const targetSubRegionIndex = meshElems.m_toElementSubRegion[iperf]; + ElementRegionBase const & region = meshLevel.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); + ElementSubRegionBase const & subRegion = region.getSubRegion< ElementSubRegionBase >( targetSubRegionIndex ); + + arrayView2d< const real64 > const perfLocation = wellSubRegion.getPerforationData()->getLocation(); + + string regionString( region.getName() ); - arrayView2d< const real64 > const perfLocation = wellSubRegion.getPerforationData()->getLocation(); + //stdVector< buffer_unit_type > rcvBuffer; + //stdVector< string > destBuffer; + //buffer_unit_type * bRegion = new buffer_unit_type[ sizeof(regionString.size()) ](); + //bufferOps::Pack< true >( bRegion, regionString ); + + stdVector< string > rcvRegionBuffer; + MpiWrapper::gather( ®ion.getName(), 1, rcvRegionBuffer.data(), 1, 0, MPI_COMM_GEOS ); + // buffer_unit_type const * localDestBuffer = rcvBuffer.data(); + // bufferOps::Unpack( localDestBuffer, destBuffer[0] ); + + if( MpiWrapper::commRank() == 0 ) + { + std::cout << "size buffer region "<< rcvRegionBuffer.size()<< " id "<< rcvRegionBuffer[0]<< std::endl; dataPerforation.addRow( iperf, wellSubRegion.getPerforationData()->getWellElements()[iperf], perfLocation[iperf], region.getName(), subRegion.getName(), cellId ); - } + } - if( dataPerforation.getCellsData().size() > 0 ) - { - TableTextFormatter const formatter( layoutPerforation ); - GEOS_LOG( formatter.toString( dataPerforation )); - } - } ); + } + if( dataPerforation.getCellsData().size() > 0 ) + { + TableTextFormatter const formatter( layoutPerforation ); + GEOS_LOG_RANK_0( formatter.toString( dataPerforation )); + } + } ); // communicate to rebuild global node info since we modified global ordering From f29df71b83982b54069abc3f0d36464e898be451 Mon Sep 17 00:00:00 2001 From: arng40 Date: Fri, 13 Feb 2026 10:45:59 +0100 Subject: [PATCH 7/8] MPI gather string - code MPI duplicate with negative pressure cell --- src/coreComponents/common/MpiWrapper.cpp | 105 ++++++++++++++++++ src/coreComponents/common/MpiWrapper.hpp | 6 + .../mesh/ElementRegionManager.cpp | 51 ++++++--- 3 files changed, 148 insertions(+), 14 deletions(-) diff --git a/src/coreComponents/common/MpiWrapper.cpp b/src/coreComponents/common/MpiWrapper.cpp index b80d0aa802a..682dbc2facc 100644 --- a/src/coreComponents/common/MpiWrapper.cpp +++ b/src/coreComponents/common/MpiWrapper.cpp @@ -518,6 +518,111 @@ template<> MPI_Datatype getMpiPairType< double, double >() } /* namespace internal */ + +template<> +void MpiWrapper::gatherStringSet< stdVector >( stdVector< string > const & strSet, + stdVector< string > & result, + MPI_Comm MPI_PARAM( comm )) +{ + int rank = MpiWrapper::commRank(); + int size = MpiWrapper::commSize(); + + string localAllStrings; + stdVector< int > localStrSizes; + localAllStrings.reserve( strSet.size() * 32 ); + + for( auto const & str : strSet ) + { + localAllStrings += str; + localStrSizes.push_back( static_cast< int >(str.size())); + } + + //1 - We gather the metadata across all ranks + struct MetaData + { + int count; + int size; + }; + MetaData localMeta = { static_cast< int >(strSet.size()), static_cast< int >(localAllStrings.size()) }; + + stdVector< MetaData > allMeta; + if( rank == 0 ) + { + allMeta.resize( size ); + } + + int localMetaArr[2] = { localMeta.count, localMeta.size }; + stdVector< int > allMetaArr; + + if( rank == 0 ) + allMetaArr.resize( size * 2 ); + + MpiWrapper::gather( localMetaArr, 2, allMetaArr.data(), 2, 0, comm ); + + //2 - Prepare displacement arrays for rank 0 + stdVector< int > allStrSizes; + string allStrings; + stdVector< int > counts_counts( size ); + stdVector< int > displs_counts( size ); + stdVector< int > counts_chars( size ); + stdVector< int > displs_chars( size ); + + int totalStrCount = 0; + + if( rank == 0 ) + { + int currentCountOffset = 0; + int currentCharOffset = 0; + + for( int currRank = 0; currRank < size; ++currRank ) + { + int c_count = allMetaArr[2*currRank]; + int c_size = allMetaArr[2*currRank + 1]; + + counts_counts[currRank] = c_count; + displs_counts[currRank] = currentCountOffset; + + counts_chars[currRank] = c_size; + displs_chars[currRank] = currentCharOffset; + + currentCountOffset += c_count; + currentCharOffset += c_size; + } + + totalStrCount = currentCountOffset; + + allStrSizes.resize( totalStrCount ); + allStrings.resize( currentCharOffset ); + } + + // 3. Gatherv des tailles de chaรฎnes + MpiWrapper::gatherv( localStrSizes.data(), localMeta.count, + allStrSizes.data(), counts_counts.data(), displs_counts.data(), + 0, comm ); + + // 4. Gatherv du contenu (chars) + MpiWrapper::gatherv( localAllStrings.c_str(), localMeta.size, + allStrings.data(), counts_chars.data(), displs_chars.data(), + 0, comm ); + + // 5. Reconstruction + if( rank == 0 ) + { + const char * dataPtr = allStrings.c_str(); + int currentOffset = 0; + for( int i = 0; i < totalStrCount; ++i ) + { + int len = allStrSizes[i]; + + std::string const temp( dataPtr + currentOffset, len ); + + result.emplace( result.end(), temp ); + currentOffset += len; + } + } + +} + } /* namespace geos */ #if defined(__clang__) diff --git a/src/coreComponents/common/MpiWrapper.hpp b/src/coreComponents/common/MpiWrapper.hpp index fd3264a5317..4b7428770d7 100644 --- a/src/coreComponents/common/MpiWrapper.hpp +++ b/src/coreComponents/common/MpiWrapper.hpp @@ -22,6 +22,7 @@ #include "common/DataTypes.hpp" #include "common/Span.hpp" +#include "common/StdContainerWrappers.hpp" #include "common/TypesHelpers.hpp" #include @@ -331,6 +332,11 @@ struct MpiWrapper */ static int nodeCommSize(); + template< template< class > class CONTAINER = stdVector > + void static gatherStringSet( CONTAINER< string > const & strSet, + CONTAINER< string > & result, + MPI_Comm MPI_PARAM( comm )); + /** * @brief Strongly typed wrapper around MPI_Allgather. * @tparam T_SEND The pointer type for \p sendbuf diff --git a/src/coreComponents/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 3011115ed4e..939782c4eab 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -14,10 +14,12 @@ */ #include +#include #include #include "ElementRegionManager.hpp" +#include "codingUtilities/RTTypes.hpp" #include "common/DataLayouts.hpp" #include "common/DataTypes.hpp" #include "common/MpiWrapper.hpp" @@ -238,6 +240,9 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); + int rank = MpiWrapper::commRank(); + // int size = MpiWrapper::commSize(); + forElementRegions< WellElementRegion >( [&]( WellElementRegion const & wellRegion ){ WellElementSubRegion const & @@ -258,28 +263,46 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM ElementRegionBase const & region = meshLevel.getElemManager().getRegion< ElementRegionBase >( targetRegionIndex ); ElementSubRegionBase const & subRegion = region.getSubRegion< ElementSubRegionBase >( targetSubRegionIndex ); - arrayView2d< const real64 > const perfLocation = wellSubRegion.getPerforationData()->getLocation(); + arrayView2d< const real64 > perfLocation = wellSubRegion.getPerforationData()->getLocation(); - string regionString( region.getName() ); + stdVector< string > setRegionName( {region.getName()} ); + stdVector< string > resultRegionName; + stdVector< string > setSubRegionName( {subRegion.getName()} ); + stdVector< string > resultSubRegionName; + MpiWrapper::gatherStringSet< stdVector >( setRegionName, resultRegionName, MPI_COMM_WORLD ); + MpiWrapper::gatherStringSet< stdVector >( setSubRegionName, resultSubRegionName, MPI_COMM_WORLD ); + array1d< globalIndex > rcvPerfo; + rcvPerfo.resize( 3 ); - //stdVector< buffer_unit_type > rcvBuffer; - //stdVector< string > destBuffer; - //buffer_unit_type * bRegion = new buffer_unit_type[ sizeof(regionString.size()) ](); - //bufferOps::Pack< true >( bRegion, regionString ); + if( rank == 0 ) + { + if( perfLocation.empty()) + { + MpiWrapper::recv( rcvPerfo, MPI_ANY_SOURCE, 123, MPI_COMM_WORLD, MPI_STATUS_IGNORE ); + } + else + { + LvArray::forValuesInSlice( perfLocation.toSlice(), [&]( real64 const & v ){rcvPerfo.emplace_back( v ); } ); + } + } + else + { + if( !perfLocation.empty()) + { + MpiWrapper::send( perfLocation.data(), 3, 0, 123, MPI_COMM_WORLD ); + } + } - stdVector< string > rcvRegionBuffer; - MpiWrapper::gather( ®ion.getName(), 1, rcvRegionBuffer.data(), 1, 0, MPI_COMM_GEOS ); - // buffer_unit_type const * localDestBuffer = rcvBuffer.data(); - // bufferOps::Unpack( localDestBuffer, destBuffer[0] ); - if( MpiWrapper::commRank() == 0 ) + if( rank == 0 ) { - std::cout << "size buffer region "<< rcvRegionBuffer.size()<< " id "<< rcvRegionBuffer[0]<< std::endl; + std::cout << "destBuffer "<< rcvPerfo[0]<< " "<getWellElements()[iperf], perfLocation[iperf], - region.getName(), - subRegion.getName(), + resultRegionName[0], + setSubRegionName[0], cellId ); } From c09f8bfe8f62cf67371ed6ce0df397496660e834 Mon Sep 17 00:00:00 2001 From: arng40 Date: Fri, 13 Feb 2026 10:49:47 +0100 Subject: [PATCH 8/8] Squashed commit of the following: MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit commit 82269bc4f0e92fc25f8c9e8ee054b6b39627597d Merge: 28643c97c6 cf245f7b94 Author: Nicola Castelletto <38361926+castelletto1@users.noreply.github.com> Date: Thu Feb 12 17:31:28 2026 -0800 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 28643c97c6ac83f73ac7eadc1cfe34e68803058e Merge: 9be8d0b3d3 f6d06635c9 Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Thu Feb 12 13:51:40 2026 +0100 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 9be8d0b3d37d32ba2b9b455493ea59b0583ffe54 Merge: fe602227ee be4591402d Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Wed Feb 11 16:52:42 2026 +0100 Merge branch 'develop' into feature/rey/negative-pressure-cells commit fe602227ee7e5f51e70ede6f28528d685ceed658 Author: MelReyCG Date: Mon Feb 9 16:14:35 2026 +0100 ๐Ÿ”Š pressures/densities equal to 0.0 are also reported (not only negatives) commit 168bd0c005e49d65a82403c580c0735fb49a0dc8 Merge: b3094e4dc7 4103efdb30 Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Tue Nov 18 11:30:20 2025 +0100 Merge branch 'develop' into feature/rey/negative-pressure-cells commit b3094e4dc703a701b11ad71eed5beb4910c1efa5 Merge: 370da84f2a 889ea4e797 Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Wed Nov 5 09:48:33 2025 +0100 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 370da84f2ac20c11617360682ca414b000a0b38e Author: MelReyCG Date: Tue Nov 4 17:10:52 2025 +0100 ๐Ÿ› unit test merge fix commit 89c5505246f650df3d3375ab054a9e7bcbe8d58c Author: MelReyCG Date: Tue Nov 4 15:17:52 2025 +0100 ๐Ÿ“ doc fix commit b5d66becd7e175c149bf8437aeec270837d08eac Author: MelReyCG Date: Tue Nov 4 14:46:27 2025 +0100 ๐ŸŽจ ๐Ÿ“ code style & docs commit b3bfe1a880ddd9690b665c862e87e363d31321e1 Author: MelReyCG Date: Tue Nov 4 14:38:43 2025 +0100 ๐ŸŽจ copilot code checks commit 69cd64c57ffc358b2706777d423ee257b86eb1e7 Author: MelReyCG Date: Wed Oct 29 10:50:29 2025 +0100 โšฐ๏ธ merge missing variable commit c46879b8dfe0cbb0ec8f1e10ad3d73340447184e Merge: 9cf5abe9c1 178874f7a8 Author: MelReyCG Date: Wed Oct 29 10:28:02 2025 +0100 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit 9cf5abe9c1f474adfb994acdf0a4fffe1589bd29 Author: MelReyCG Date: Tue Oct 28 14:36:31 2025 +0100 ๐Ÿ› typo 2 commit b78a3d79c1cfab6cbd828942f48f0144106d05a7 Author: MelReyCG Date: Mon Oct 27 16:52:03 2025 +0100 ๐Ÿ“ small precision on loginfo commit 79e1b11134179b39bb9fcefe0ac520423b379ddf Author: MelReyCG Date: Tue Oct 28 14:22:50 2025 +0100 ๐Ÿ› typo commit 2a45e9867fd9ec1acf659c5fcb799138d6dd8c57 Merge: e4431512e6 7416c57d6b Author: MelReyCG Date: Tue Oct 28 11:54:26 2025 +0100 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit e4431512e691e9a748d289939921466adc934fd9 Author: MelReyCG Date: Tue Aug 26 15:46:38 2025 +0200 โ™ป๏ธ renamings (after Pavel review) commit 7812f759891d220ae496cb3e8fb193b56a52f065 Merge: f18cb8ec26 0cb3f1fec7 Author: MelReyCG Date: Tue Aug 26 10:46:30 2025 +0200 Merge branch 'develop' into feature/rey/negative-pressure-cells commit f18cb8ec26f4d5d3f8f0aff69d14f73ecdf9dcd7 Merge: de56d0df5f 43216af935 Author: Pavel Tomin Date: Thu Aug 7 12:43:10 2025 -0500 Merge branch 'develop' into feature/rey/negative-pressure-cells commit de56d0df5ff48a2037ddc3424200be8b2ac9592d Author: Pavel Tomin Date: Wed Aug 6 17:16:17 2025 -0500 code style commit a32f7d700dbcc4b2a12df9ee0d75ed4224e8f5fc Merge: 8bf347068d b7e609f00a Author: Pavel Tomin Date: Wed Aug 6 11:54:24 2025 -0500 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 8bf347068d2c1ac95b60cd5d4c3429e1ba7eeec5 Author: MelReyCG Date: Mon Jul 21 16:01:43 2025 +0200 ๐Ÿ› compil fix commit 96ffe23a6e2ce8e2025dc785337f2e7f77b73aa6 Merge: 7a83961353 b05541b588 Author: MelReyCG Date: Mon Jul 21 16:01:31 2025 +0200 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit 7a839613533aaa907f0d2354350a792761aae8f7 Author: MelReyCG Date: Fri Jul 18 15:42:12 2025 +0200 ๐Ÿ“ added an idea commit 2e4193cef7744bca39f34872d845e5389f1cd30c Author: MelReyCG Date: Fri Jul 18 15:33:57 2025 +0200 ๐Ÿ“ previous commit fix commit 409f4619c648d9bd27d38c3ad10671e3c3574b34 Author: MelReyCG Date: Fri Jul 18 14:24:51 2025 +0200 ๐Ÿ“ documentation updates commit 937c446b42e21751566dc2c4076982cbc433c55f Author: MelReyCG Date: Thu Jul 17 16:12:09 2025 +0200 โšฐ๏ธ removed dead code commit aeb8ea02052759bcd5fee2a6a0e7088e5f6201f9 Author: MelReyCG Date: Wed Jul 16 17:21:42 2025 +0200 ๐Ÿ“ฆ schema commit 4b46f19f10e6e87825174dd44d1050c2240d7f37 Author: MelReyCG Date: Wed Jul 16 17:21:18 2025 +0200 ๐Ÿ› solve a bug where the last line of the table was cut (when the r0 had no content) commit 92e9abbdb5be474697c206a8c7e524b20a15b60e Author: MelReyCG Date: Wed Jul 16 17:04:55 2025 +0200 ๐Ÿงช adding a (failing) test to highlight a bug commit 63fe7e79c3188af554f8ac48aa6a9c0d23317e49 Merge: 8d50b4cb02 aa33980b27 Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Wed Jul 16 16:13:07 2025 +0200 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 8d50b4cb020fd9bd16c75289abbe22582df7e41e Author: MelReyCG Date: Wed Jul 16 16:11:18 2025 +0200 โ™ป๏ธ set constant params const commit 89fe09444f9eed933eecd6b7634c94d3a628304a Merge: ef29d771be 7e85e747ec Author: MelReyCG Date: Tue Jul 8 15:50:59 2025 +0200 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit ef29d771be9d734f01483d74530f682165c55fc6 Author: MelReyCG Date: Fri Jul 4 15:55:31 2025 +0200 ranksStrsDisps -> ranksStrsOffsets commit 61c4199041e44e81be0198b5d5e84e1db4478b4a Author: MelReyCG Date: Fri Jul 4 15:38:45 2025 +0200 ๐Ÿ“ doc fix commit 5af2601ebf61f44c61ce5c6146e85e9ca2c14641 Author: MelReyCG Date: Wed Jul 2 12:05:25 2025 +0200 ๐Ÿ“ missing last docs commit de75861c879c5c2a6d9d7272c16abba3128853d9 Merge: 8eb526bf24 125ffb92a8 Author: MelReyCG <122801580+MelReyCG@users.noreply.github.com> Date: Wed Jul 2 11:52:59 2025 +0200 Merge branch 'develop' into feature/rey/negative-pressure-cells commit 8eb526bf24fc2c9908cf35281fb70cd859f9c9a7 Author: MelReyCG Date: Wed Jul 2 11:51:48 2025 +0200 ๐ŸŽจ uncrustify commit 411636da54c1ab1acde6cd359b9e7bea006df4d1 Author: MelReyCG Date: Wed Jul 2 11:48:22 2025 +0200 ๐Ÿ“ Adding las docs commit e5ae68ec78694ae737a23292ed8581451b149f4d Author: MelReyCG Date: Tue Jul 1 17:54:36 2025 +0200 โœ… adding mpi tables unit test commit 112be17fb935f425d1461fbd1ceb99e4ab292e10 Author: MelReyCG Date: Tue Jul 1 17:54:09 2025 +0200 ๐Ÿ“ updating documentation commit b086a4f90e2e275edc7fccae615da50d2e98cf02 Author: MelReyCG Date: Tue Jul 1 15:46:33 2025 +0200 ๐Ÿ› fixing scarce crash when mpi-tables are constructed from more than 2 ranks commit 04754361d56ad26a96486c5bb896b33c3a3b8a54 Author: MelReyCG Date: Tue Jul 1 15:45:59 2025 +0200 โšฐ๏ธ unused variables commit 72e8ecdf301da46d80354a3e045289e8a5e9917d Author: MelReyCG Date: Fri Jun 6 10:22:27 2025 +0200 ๐Ÿ’„ transposing table layout for clarity (user review) commit 390e451f746fb6d3612b93955560d5f209f58031 Author: MelReyCG Date: Fri Jun 6 10:18:17 2025 +0200 โ™ป๏ธ Adding missing signatures commit 9b6cd4d79a3ea03ee1c28507303f54cf432ac053 Merge: 5ab8421400 c2768e6ecd Author: MelReyCG Date: Wed Jun 4 11:46:45 2025 +0200 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit 5ab8421400cf2143049706dc84eb022be131f910 Author: MelReyCG Date: Wed Jun 4 11:46:15 2025 +0200 โœจ adding ranks separator titles commit 6a43cb9dd0d24def34d0acad25ed376cf1c0a8df Author: MelReyCG Date: Tue Jun 3 18:24:38 2025 +0200 โœจ finishing MPI tables with a different approach (log output on rank0 only) commit 6487d6c2b3e774eb86e0ea64f4a6f389cb99e6f2 Author: MelReyCG Date: Mon Jun 2 14:13:33 2025 +0200 โœจ ๐Ÿงช first attempt at creating MPI tables commit 854dc3e11a464e60b365edf729b1c09a06172c64 Author: MelReyCG Date: Wed May 28 14:34:17 2025 +0200 โ™ป๏ธ minor refactor of TableFormatter.cpp commit 66f2d727e9c1df6c1c5c805082e5be0f76db5657 Author: MelReyCG Date: Wed May 28 14:33:57 2025 +0200 โœจ adding info to warn the user to increase a logLevel to get the report commit 4a9f7803e20fb9cb676864b68eea3f5863bbea52 Author: MelReyCG Date: Wed May 28 14:16:22 2025 +0200 โœจextending reported data to pressure/density (IdReporter -> ElementReporter) commit 00d8c577c219a1785c1ebe581838d81a08719119 Author: MelReyCG Date: Wed May 28 11:00:31 2025 +0200 โšฐ๏ธ dead code commit da89daf1784259f994b7be243ff848d54089af4b Author: MelReyCG Date: Tue May 27 17:09:02 2025 +0200 ๐Ÿ’„ new table format commit e7f73cf188698b44b110f98ef235314c50d5e4c0 Author: MelReyCG Date: Mon May 26 17:49:19 2025 +0200 โœจ offering control on alignement when using columns-free table layouts commit 3fd0c12e336f924a27a00b7d273b81fdab774c67 Author: MelReyCG Date: Mon May 26 17:48:52 2025 +0200 โ™ป๏ธ ๐Ÿ› simplifying & fixing visible columns counting commit 8f68eb515f0534be35d467dfc40a3a42b52c3e06 Author: MelReyCG Date: Mon May 26 17:47:25 2025 +0200 โ™ป๏ธ refactoring TableFormatter: give control to inheriting classes commit 5744a2af42292f5ed05e6cd88f289150019a87dd Merge: 58d1dce874 ef5e9406d7 Author: MelReyCG Date: Mon May 26 14:06:10 2025 +0200 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit 58d1dce87471a359eb18a982f548e4c0dfb0812a Author: MelReyCG Date: Mon May 26 11:52:38 2025 +0200 โ™ป๏ธ refactoring TableFormatter: give control to inheriting classes commit 3eb2293a826a844df7c0a7f8b678f165fd0f6c79 Author: MelReyCG Date: Mon May 26 11:48:38 2025 +0200 โ™ป๏ธ removing double assessor + exposing non-const version commit 5fc7ea63c7da07b4f7364e9d2a5035c052a858ef Author: MelReyCG Date: Fri May 23 12:28:09 2025 +0200 โœจ adding support for no column titled table layout commit 47f5c5ce91bd396d355e4034b10caa26428ab74d Author: MelReyCG Date: Fri May 23 11:09:44 2025 +0200 โœจ adding table indentation commit de506c60bb0ee7b09eaac5518f72efe2a730dbb8 Author: MelReyCG Date: Fri May 23 10:09:42 2025 +0200 โœจ adding table formatting to allow for showing more data commit 0803a92e910f5274258ad742fd227cec4e5b9151 Author: MelReyCG Date: Thu May 22 16:53:04 2025 +0200 ๐Ÿ› CUDA crash fix + bug fix (data not moved from device) commit d29db451bc30b9c76a90873b3b0910de135d859e Author: MelReyCG Date: Tue May 20 17:42:52 2025 +0200 ๐Ÿ› adding one more barrier to not get the msg cut by other msgs. commit e0b104d7c12ac6230138e9d04b070684c87f5ad2 Author: MelReyCG Date: Tue May 20 17:22:31 2025 +0200 ๐Ÿ’„ msg slight rewriting commit 3f97e2bd5527ba9269893344374492be5982871b Author: MelReyCG Date: Tue May 20 14:33:40 2025 +0200 ๐Ÿ› fix for a CUDA target: explicit constructor call commit 6766e7e6a682bf49d6e99ba9a91b63ed0b6be961 Merge: 056ac9cf89 5a65081dca Author: MelReyCG Date: Tue May 20 11:37:03 2025 +0200 Merge remote-tracking branch 'origin/develop' into feature/rey/negative-pressure-cells commit 056ac9cf89e5d9d918a870c7525d734e896e9288 Author: MelReyCG Date: Tue May 20 11:29:21 2025 +0200 โœจ implementation of neg pressure ids output on other models commit 18703e3ef100a37c09ea9b8a057c9949a3e4b25a Author: MelReyCG Date: Tue May 20 11:28:20 2025 +0200 โ™ป๏ธ refactors to prevent computing mpi reduction twice + constness + naming commit 14a546fbb1eb3ef5254207058be2c63fbb2a8fb6 Author: MelReyCG Date: Mon May 19 16:10:47 2025 +0200 ๐Ÿ’„ adding msg precisions commit 4ebe3308a3e1486c85411c2ecfaf7ae981fb24e0 Author: MelReyCG Date: Mon May 19 15:48:26 2025 +0200 ๐Ÿ’„ slight msg change / fix 2 commit 5cf681482329ec2fd28aad90acc2a5d5b1160d0a Author: MelReyCG Date: Mon May 19 15:40:18 2025 +0200 ๐Ÿ’„ slight msg change / fix commit 96cd0a844cabf0e7465fc9d046edee117890e9ae Author: MelReyCG Date: Fri May 16 18:10:08 2025 +0200 โœจ activation de l'output des valeurs problรฉmatique en fยฐ d'un nouveau logLevel "SolutionDetails" commit cb0ae03efeef7fff41d2378e11e97d1c254174ca Author: MelReyCG Date: Fri May 16 12:34:53 2025 +0200 โ™ป๏ธ proper cpp/hpp file repartition commit 1434bd389f2db1f2cc9b060e798baae5a0cd9cfb Author: MelReyCG Date: Fri May 16 12:09:57 2025 +0200 โ™ป๏ธ removing unnecessary template commit 8d158fc0ded2c0dbf182e1317321a6aeeace65e5 Author: MelReyCG Date: Fri May 16 11:56:59 2025 +0200 โœจ added the IdReporter classes to facilitate outputing ids & their count from kernels (optionally) commit 76f4644c9011231612342e0d3e38575bec9891e8 Author: MelReyCG Date: Wed May 7 18:24:23 2025 +0200 ๐Ÿ’ฉ first attempt at outputing wrong cells in multiphase cases... work on CPU but doesn't on GPU build + suboptimal --- src/coreComponents/common/CMakeLists.txt | 10 +- .../common/format/table/TableData.cpp | 10 - .../common/format/table/TableData.hpp | 21 +- .../common/format/table/TableFormatter.cpp | 70 ++++--- .../common/format/table/TableFormatter.hpp | 62 ++++-- .../common/format/table/TableLayout.cpp | 43 +++- .../common/format/table/TableLayout.hpp | 85 ++++++-- .../format/table/TableMpiComponents.cpp | 160 +++++++++++++++ .../format/table/TableMpiComponents.hpp | 110 ++++++++++ .../format/table/unitTests/CMakeLists.txt | 20 ++ .../format/table/unitTests/testMpiTable.cpp | 149 ++++++++++++++ .../format/table/unitTests/testTable.cpp | 26 ++- .../parameters/KValueFlashParameters.cpp | 4 +- .../physicsSolvers/LogLevelsInfo.hpp | 6 + .../physicsSolvers/fluidFlow/CMakeLists.txt | 3 + .../fluidFlow/CompositionalMultiphaseFVM.cpp | 117 ++++++----- .../CompositionalMultiphaseHybridFVM.cpp | 4 +- .../fluidFlow/SinglePhaseBase.cpp | 32 +-- .../fluidFlow/SolutionCheckHelpers.cpp | 118 +++++++++++ .../fluidFlow/SolutionCheckHelpers.hpp | 190 ++++++++++++++++++ .../kernels/SolutionCheckKernelsHelpers.hpp | 137 +++++++++++++ .../compositional/SolutionCheckKernel.hpp | 159 +++++++++------ .../ThermalSolutionCheckKernel.hpp | 12 +- .../singlePhase/SolutionCheckKernel.hpp | 33 +-- .../wells/CompositionalMultiphaseWell.cpp | 118 +++++------ .../fluidFlow/wells/SinglePhaseWell.cpp | 35 ++-- 26 files changed, 1415 insertions(+), 319 deletions(-) create mode 100644 src/coreComponents/common/format/table/TableMpiComponents.cpp create mode 100644 src/coreComponents/common/format/table/TableMpiComponents.hpp create mode 100644 src/coreComponents/common/format/table/unitTests/testMpiTable.cpp create mode 100644 src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.cpp create mode 100644 src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp create mode 100644 src/coreComponents/physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp diff --git a/src/coreComponents/common/CMakeLists.txt b/src/coreComponents/common/CMakeLists.txt index 0e37ba56703..d0143b9564b 100644 --- a/src/coreComponents/common/CMakeLists.txt +++ b/src/coreComponents/common/CMakeLists.txt @@ -24,9 +24,10 @@ Also provides commonly used components for such as logging, formatting, memory a # set( common_headers ${CMAKE_BINARY_DIR}/include/common/GeosxConfig.hpp - format/table/TableLayout.hpp - format/table/TableFormatter.hpp format/table/TableData.hpp + format/table/TableFormatter.hpp + format/table/TableLayout.hpp + format/table/TableMpiComponents.hpp format/EnumStrings.hpp format/LogPart.hpp format/Format.hpp @@ -71,9 +72,10 @@ endif( ) # Specify all sources # set( common_sources - format/table/TableLayout.cpp - format/table/TableFormatter.cpp format/table/TableData.cpp + format/table/TableFormatter.cpp + format/table/TableLayout.cpp + format/table/TableMpiComponents.cpp format/LogPart.cpp format/StringUtilities.cpp logger/GeosExceptions.cpp diff --git a/src/coreComponents/common/format/table/TableData.cpp b/src/coreComponents/common/format/table/TableData.cpp index 072e3b1ca79..b9bfc885584 100644 --- a/src/coreComponents/common/format/table/TableData.cpp +++ b/src/coreComponents/common/format/table/TableData.cpp @@ -89,16 +89,6 @@ void TableData::clear() getErrorsList().clear(); } -stdVector< stdVector< TableData::CellData > > const & TableData::getTableDataRows() const -{ - return m_rows; -} - -stdVector< stdVector< TableData::CellData > > & TableData::getTableDataRows() -{ - return m_rows; -} - void TableData2D::collectTableValues( arrayView1d< real64 const > dim0AxisCoordinates, arrayView1d< real64 const > dim1AxisCoordinates, arrayView1d< real64 const > values, diff --git a/src/coreComponents/common/format/table/TableData.hpp b/src/coreComponents/common/format/table/TableData.hpp index 22135105ecf..e36afb037c5 100644 --- a/src/coreComponents/common/format/table/TableData.hpp +++ b/src/coreComponents/common/format/table/TableData.hpp @@ -111,16 +111,6 @@ class TableData void clearErrors() { m_errors->clear(); } - /** - * @return The const rows of the table - */ - stdVector< stdVector< CellData > > const & getTableDataRows() const; - - /** - * @return The rows of the table - */ - stdVector< stdVector< CellData > > & getTableDataRows(); - /** * @brief Get all error messages * @return The vector of error messages @@ -133,16 +123,19 @@ class TableData DataRows const & getCellsData() const { return m_rows; } + /** + * @return The const table data rows + */ + DataRows & getCellsData() + { return m_rows; } + /** * @brief Comparison operator for data rows * @param comparingTable The tableData values to compare * @return The comparison result */ inline bool operator==( TableData const & comparingTable ) const - { - - return getCellsData() == comparingTable.getCellsData(); - } + { return getCellsData() == comparingTable.getCellsData(); } /** * @brief Get all error messages diff --git a/src/coreComponents/common/format/table/TableFormatter.cpp b/src/coreComponents/common/format/table/TableFormatter.cpp index e34eaa8d37c..bf5779fdc23 100644 --- a/src/coreComponents/common/format/table/TableFormatter.cpp +++ b/src/coreComponents/common/format/table/TableFormatter.cpp @@ -18,10 +18,10 @@ * @file TableFormatter.cpp */ -#include "TableFormatter.hpp" #include #include "common/format/StringUtilities.hpp" #include "common/logger/Logger.hpp" +#include "TableFormatter.hpp" namespace geos { @@ -157,7 +157,8 @@ string TableCSVFormatter::headerToString() const string TableCSVFormatter::dataToString( TableData const & tableData ) const { - RowsCellInput const rowsValues( tableData.getTableDataRows() ); + + RowsCellInput const rowsValues( tableData.getCellsData() ); string result; size_t total_size = 0; for( auto const & row : rowsValues ) @@ -232,12 +233,13 @@ string TableTextFormatter::toString< TableData >( TableData const & tableData ) initalizeTableGrids( m_tableLayout, tableData, headerCellsLayout, dataCellsLayout, errorCellsLayout, - tableTotalWidth ); - outputTable( m_tableLayout, tableOutput, - headerCellsLayout, dataCellsLayout, errorCellsLayout, - tableTotalWidth ); + tableTotalWidth, nullptr ); + + string const sepLine = string( tableTotalWidth, m_horizontalLine ); + outputTableHeader( tableOutput, m_tableLayout, headerCellsLayout, sepLine ); + outputTableData( tableOutput, m_tableLayout, dataCellsLayout ); + outputTableFooter( tableOutput, m_tableLayout, errorCellsLayout, sepLine, !dataCellsLayout.empty() ); - getErrorsList().clear(); return tableOutput.str(); } @@ -246,10 +248,11 @@ void TableTextFormatter::initalizeTableGrids( PreparedTableLayout const & tableL CellLayoutRows & headerCellsLayout, CellLayoutRows & dataCellsLayout, CellLayoutRows & errorCellsLayout, - size_t & tableTotalWidth ) const + size_t & tableTotalWidth, + ColumnWidthModifier columnWidthModifier ) const { + RowsCellInput const & inputDataValues( tableInputData.getCellsData() ); bool const hasColumnLayout = tableLayout.getColumnLayersCount() > 0; - RowsCellInput const & inputDataValues( tableInputData.getTableDataRows() ); size_t const inputDataRowsCount = !inputDataValues.empty() ? inputDataValues.front().size() : 0; size_t const nbVisibleColumns = std::max( size_t( 1 ), ( hasColumnLayout ? tableLayout.getVisibleLowermostColumnCount() : @@ -282,6 +285,9 @@ void TableTextFormatter::initalizeTableGrids( PreparedTableLayout const & tableL stretchColumnsByMergedCellsWidth( columnsWidth, dataCellsLayout, tableLayout, true ); stretchColumnsByMergedCellsWidth( columnsWidth, errorCellsLayout, tableLayout, true ); + if( columnWidthModifier ) + columnWidthModifier( columnsWidth ); + // the columns width array is now sized after all the table, we can compute the total table width tableTotalWidth = tableLayout.getBorderMargin() * 2 + 2; for( size_t columnId = 0; columnId < columnsWidth.size(); ++columnId ) @@ -308,14 +314,14 @@ void TableTextFormatter::populateTitleCellsLayout( PreparedTableLayout const & t // the title row consists in a row of cells merging with the last cell containing the title text headerCellsLayout.emplace_back() = { stdVector< TableLayout::CellLayout >( nbVisibleColumns, - TableLayout::CellLayout( CellType::MergeNext ) ), // cells + TableLayout::CellLayout( CellType::MergeNext ) ), // cells titleInput.getHeight(), // sublinesCount }; headerCellsLayout.back().cells.back() = titleInput; headerCellsLayout.emplace_back() = { stdVector< TableLayout::CellLayout >( nbVisibleColumns, - TableLayout::CellLayout( CellType::Separator ) ), // cells + TableLayout::CellLayout( CellType::Separator ) ), // cells 1, // sublinesCount }; } @@ -494,8 +500,8 @@ void TableTextFormatter::populateDataCellsLayout( PreparedTableLayout const & ta string_view( &m_horizontalLine, 1 ) : string_view( inputCell.value ); TableLayout::Alignment const alignment = inputCell.type == CellType::Header ? - tableLayout.defaultHeaderAlignment : - tableLayout.defaultValueAlignment; + tableLayout.getDefaultHeaderAlignment() : + tableLayout.getDefaultValueAlignment(); TableLayout::CellLayout & outputCell = outputRow.cells[idxColumn]; outputCell = TableLayout::CellLayout( inputCell.type, alignment ); @@ -702,34 +708,44 @@ void TableTextFormatter::applyColumnsWidth( stdVector< size_t > const & columnsW } } -void TableTextFormatter::outputTable( PreparedTableLayout const & tableLayout, - std::ostream & tableOutput, - CellLayoutRows const & headerCellsLayout, - CellLayoutRows const & dataCellsLayout, - CellLayoutRows & errorCellsLayout, - size_t const tableTotalWidth ) const +void TableTextFormatter::outputTableHeader( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & headerCellsLayout, + string_view sepLine ) const { - string const sepLine = string( tableTotalWidth, m_horizontalLine ); if( tableLayout.isLineBreakEnabled()) { tableOutput << '\n'; } - tableOutput << sepLine << '\n'; + tableOutput << tableLayout.getIndentationStr() << sepLine << '\n'; outputLines( tableLayout, headerCellsLayout, tableOutput ); +} - if( !dataCellsLayout.empty()) +void TableTextFormatter::outputTableData( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & dataCellsLayout ) const +{ + if( !dataCellsLayout.empty() ) { outputLines( tableLayout, dataCellsLayout, tableOutput ); } +} +void TableTextFormatter::outputTableFooter( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows & errorCellsLayout, + string_view sepLine, + bool hasData ) const +{ if( !errorCellsLayout.empty()) { outputErrors( tableLayout, errorCellsLayout, tableOutput ); } - if( !dataCellsLayout.empty() || getErrorsList().hasErrors()) - tableOutput << sepLine; - + if( hasData || !errorCellsLayout.empty() ) + { + tableOutput << tableLayout.getIndentationStr() << sepLine; + } if( tableLayout.isLineBreakEnabled()) { @@ -816,6 +832,7 @@ void TableTextFormatter::outputLines( PreparedTableLayout const & tableLayout, if( isLeftBorderCell ) { // left table border isLeftBorderCell=false; + tableOutput << tableLayout.getIndentationStr(); tableOutput << m_verticalLine << string( nbBorderSpaces, cellSpaceChar ); } else @@ -845,4 +862,5 @@ void TableTextFormatter::outputLines( PreparedTableLayout const & tableLayout, idxRow++; } } -} + +} /* namespace geos */ diff --git a/src/coreComponents/common/format/table/TableFormatter.hpp b/src/coreComponents/common/format/table/TableFormatter.hpp index c2035bc2a1c..eed907975c2 100644 --- a/src/coreComponents/common/format/table/TableFormatter.hpp +++ b/src/coreComponents/common/format/table/TableFormatter.hpp @@ -86,7 +86,7 @@ class TableFormatter }; /** - * @brief class for CSV formatting + * @brief Class to format data in a formatted CSV format */ class TableCSVFormatter final : public TableFormatter { @@ -189,9 +189,10 @@ string TableCSVFormatter::toString< TableData >( TableData const & tableData ) c /** - * @brief class for log formatting + * @brief Class to format data in a formatted text format + * (for log output typically, expecting fixed character size). */ -class TableTextFormatter final : public TableFormatter +class TableTextFormatter : public TableFormatter { public: @@ -243,13 +244,15 @@ class TableTextFormatter final : public TableFormatter void toStream( std::ostream & outputStream, DATASOURCE const & tableData ) const { toStreamImpl( outputStream, toString( tableData ) ); } -private: +protected: /// symbol for separator construction static constexpr char m_verticalLine = '|'; /// for the extremity of a row static constexpr char m_horizontalLine = '-'; + /// A functor which allow to customize the columns width after their computation. + using ColumnWidthModifier = std::function< void ( stdVector< size_t > & ) >; /** * @brief Initializes the table layout with the given table data and prepares necessary layouts for headers and data cells. @@ -258,30 +261,54 @@ class TableTextFormatter final : public TableFormatter * @param headerCellsLayout A reference to a `CellLayoutRows` where the header cells will be populated. * @param dataCellsLayout A reference to a `CellLayoutRows` where the data cells will be populated. * @param errorCellsLayout A reference to a `CellLayoutRows` where the error cells will be populated. - * @param separatorLine A string that will be used as the table separator line + * @param tableTotalWidth A string that will be used as the table separator line + * @param columnWidthModifier A functor which allow to customize the columns width after their computation. */ void initalizeTableGrids( PreparedTableLayout const & tableLayout, TableData const & tableData, CellLayoutRows & dataCellsLayout, CellLayoutRows & headerCellsLayout, CellLayoutRows & errorCellsLayout, - size_t & tableTotalWidth ) const; + size_t & tableTotalWidth, + ColumnWidthModifier columnWidthModifier ) const; + + /** + * @brief Outputs the top part of the formatted table to the provided output stream. + * @param tableOutput A reference to an `std::ostream` where the formatted table will be written. + * @param tableLayout The layout of the table + * @param headerCellsLayout The header rows in a grid layout + * @param separatorLine A string that will be used as the table separator line + */ + void outputTableHeader( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & headerCellsLayout, + string_view separatorLine ) const; /** - * @brief Outputs the formatted table to the provided output stream. + * @brief Outputs the data part of the formatted table to the provided output stream. + * @param tableOutput A reference to an `std::ostream` where the formatted table will be written. * @param tableLayout The layout of the table + * @param dataCellsLayout The data rows in a grid layout + */ + void outputTableData( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & dataCellsLayout ) const; + + /** + * @brief Outputs the bottom part of the formatted table to the provided output stream. * @param tableOutput A reference to an `std::ostream` where the formatted table will be written. - * @param headerCellsLayout The layout of the header rows - * @param dataCellsLayout The layout of the data rows + * @param tableLayout The layout of the table + * @param separatorLine A string that will be used as the table separator line * @param errorCellsLayout The layout of the error rows - * @param separatorLine The string to be used as the table separator line + * @param hasData Indicates whether there is data in the table TableData. */ - void outputTable( PreparedTableLayout const & tableLayout, - std::ostream & tableOutput, - CellLayoutRows const & headerCellsLayout, - CellLayoutRows const & dataCellsLayout, - CellLayoutRows & errorCellsLayout, - size_t tableTotalWidth ) const; + void outputTableFooter( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows & errorCellsLayout, + string_view separatorLine, + bool hasData ) const; + +private: /** * @brief Outputs the formatted table lines to the output stream. @@ -312,7 +339,7 @@ class TableTextFormatter final : public TableFormatter */ void populateTitleCellsLayout( PreparedTableLayout const & tableLayout, CellLayoutRows & headerCellsLayout, - size_t const nbVisibleColumn ) const; + size_t nbVisibleColumn ) const; /** * @brief Populate a grid of CellLayout with all visible columns of the given table layout. @@ -402,6 +429,7 @@ class TableTextFormatter final : public TableFormatter void formatCell( std::ostream & tableOutput, TableLayout::CellLayout const & cell, size_t idxLine ) const; + }; /** diff --git a/src/coreComponents/common/format/table/TableLayout.cpp b/src/coreComponents/common/format/table/TableLayout.cpp index 55446491cd2..67401ff9e67 100644 --- a/src/coreComponents/common/format/table/TableLayout.cpp +++ b/src/coreComponents/common/format/table/TableLayout.cpp @@ -24,31 +24,41 @@ namespace geos { -void TableLayout::addColumns( stdVector< string > const & columnNames ) +TableLayout & TableLayout::addColumns( stdVector< string > const & columnNames ) { for( auto const & columnName : columnNames ) { addColumn( columnName ); } + return *this; } -void TableLayout::addColumns( stdVector< TableLayout::Column > const & columns ) +TableLayout & TableLayout::addColumns( stdVector< TableLayout::Column > const & columns ) { for( auto const & column : columns ) { addColumn( column ); } + return *this; +} + +TableLayout & TableLayout::addColumns( TableLayoutArgs columns ) +{ + processArguments( columns ); + return *this; } -void TableLayout::addColumn( string_view columnName ) +TableLayout & TableLayout::addColumn( string_view columnName ) { TableLayout::Column column = TableLayout::Column().setName( columnName ); m_tableColumns.emplace_back( column ); + return *this; } -void TableLayout::addColumn( TableLayout::Column const & column ) +TableLayout & TableLayout::addColumn( TableLayout::Column const & column ) { m_tableColumns.emplace_back( column ); + return *this; } TableLayout & TableLayout::setTitle( string_view title ) @@ -78,6 +88,23 @@ TableLayout & TableLayout::setMaxColumnWidth( size_t width ) return *this; } +TableLayout & TableLayout::setIndentation( size_t spacesCount ) +{ + m_indentation = spacesCount; + return *this; +} + +TableLayout & TableLayout::setDefaultHeaderAlignment( TableLayout::Alignment alignment ) +{ + m_defaultHeaderAlignment = alignment; + return *this; +} +TableLayout & TableLayout::setDefaultValueAlignment( TableLayout::Alignment alignment ) +{ + m_defaultValueAlignment = alignment; + return *this; +} + bool TableLayout::isLineBreakEnabled() const { return m_lineBreakAtBegin; } @@ -159,7 +186,7 @@ void TableLayout::Cell::setText( string_view text ) } TableLayout::Column::Column(): - m_header( CellType::Header, defaultHeaderAlignment ) + m_header( CellType::Header, Alignment::center ) {} TableLayout::Column::Column( string_view name, TableLayout::ColumnAlignement alignment ): @@ -298,14 +325,16 @@ TableLayout::DeepFirstIterator TableLayout::beginDeepFirst() const PreparedTableLayout::PreparedTableLayout( ): TableLayout(), m_columnLayersCount( 0 ), - m_visibleLowermostColumnCount( 0 ) + m_visibleLowermostColumnCount( 0 ), + m_indentationStr( m_indentation, ' ' ) {} PreparedTableLayout::PreparedTableLayout( TableLayout const & other ): TableLayout( other ), m_columnLayersCount( 0 ), m_totalLowermostColumnCount( 0 ), - m_visibleLowermostColumnCount( 0 ) + m_visibleLowermostColumnCount( 0 ), + m_indentationStr( m_indentation, ' ' ) { prepareLayoutRecusive( m_tableColumns, 0 ); diff --git a/src/coreComponents/common/format/table/TableLayout.hpp b/src/coreComponents/common/format/table/TableLayout.hpp index 4ef156ad760..e08474aa064 100644 --- a/src/coreComponents/common/format/table/TableLayout.hpp +++ b/src/coreComponents/common/format/table/TableLayout.hpp @@ -41,12 +41,6 @@ class TableLayout /// Type of aligment for a column enum Alignment { right, left, center }; - /// default value for columns header cells alignement - static constexpr Alignment defaultHeaderAlignment = Alignment::center; - - /// default value for data cells alignement - static constexpr Alignment defaultValueAlignment = Alignment::right; - /// Space to apply between all data and border enum MarginValue : integer { @@ -67,9 +61,9 @@ class TableLayout struct ColumnAlignement { /// Alignment for column name. By default aligned to center - Alignment headerAlignment = defaultHeaderAlignment; + Alignment headerAlignment = Alignment::center; /// Alignment for column values. By default aligned to right side - Alignment valueAlignment = defaultValueAlignment; + Alignment valueAlignment = Alignment::right; }; /** @@ -618,6 +612,18 @@ class TableLayout string_view getTitleStr() const { return m_tableTitleStr; } + /** + * @return the default value for columns header cells alignement. Used with column-free layout. + */ + Alignment getDefaultHeaderAlignment() const + { return m_defaultHeaderAlignment; } + + /** + * @return the default value for data cells alignement. Used with column-free layout. + */ + Alignment getDefaultValueAlignment() const + { return m_defaultValueAlignment; } + /** * @param title The table title * @return The tableLayout reference @@ -645,6 +651,27 @@ class TableLayout */ TableLayout & setMaxColumnWidth( size_t width ); + /** + * @brief Set the indentation of the whole table. + * @param spacesCount The number of indentation spaces. + * @return the TableLayout instance, for builder pattern + */ + TableLayout & setIndentation( size_t spacesCount ); + + /** + * @brief Sets the default value for columns header cells alignement. Used with column-free layout. + * @param alignment The desired alignment + * @return the TableLayout instance, for builder pattern + */ + TableLayout & setDefaultHeaderAlignment( Alignment alignment ); + + /** + * @brief Sets the default value for data cells alignement. Used with column-free layout. + * @param alignment The desired alignment + * @return the TableLayout instance, for builder pattern + */ + TableLayout & setDefaultValueAlignment( Alignment alignment ); + /** * @brief check if a column max width has been set * @return Truef a column max width has been set, otherwise false @@ -681,29 +708,46 @@ class TableLayout size_t const & getMaxColumnWidth() const { return m_maxColumnWidth; } + /** + * @return The number of spaces at the left of the table. + */ + size_t const & getIndentation() const + { return m_indentation; } + /** * @brief Create and add columns to the columns vector given a string vector * @param columnNames The columns name + * @return the TableLayout instance, for builder pattern */ - void addColumns( stdVector< TableLayout::Column > const & columnNames ); + TableLayout & addColumns( stdVector< Column > const & columnNames ); /** * @brief Create and add columns to the columns vector given a string vector * @param columns The columns list + * @return the TableLayout instance, for builder pattern + */ + TableLayout & addColumns( stdVector< string > const & columns ); + + /** + * @brief Create and add columns to the columns vector given a string and/or columns + * @param columns brace enclosed parameters, consisting of column names or Column instances + * @return the TableLayout instance, for builder pattern */ - void addColumns( stdVector< string > const & columns ); + TableLayout & addColumns( TableLayoutArgs columns ); /** * @brief Create and add a column to the columns vector given a string * @param columnName The column name + * @return the TableLayout instance, for builder pattern */ - void addColumn( string_view columnName ); + TableLayout & addColumn( string_view columnName ); /** * @brief Create and add a column to the columns vector given a Column * @param column Vector containing addition information on the column + * @return the TableLayout instance, for builder pattern */ - void addColumn( TableLayout::Column const & column ); + TableLayout & addColumn( Column const & column ); protected: @@ -756,6 +800,15 @@ class TableLayout /// The number of margin spaces around contents. integer m_marginValue; + /// The number of spaces at the left of the table. + size_t m_indentation = 0; + + /// default value for columns header cells alignement. + Alignment m_defaultHeaderAlignment = Alignment::center; + + /// default value for data cells alignement. + Alignment m_defaultValueAlignment = Alignment::right; + }; /** @@ -814,6 +867,12 @@ class PreparedTableLayout : public TableLayout size_t getTotalLowermostColumnCount() const { return m_totalLowermostColumnCount; } + /** + * @return A string with the correct indentation space count to precede each lines of the formatted table. + */ + string_view getIndentationStr() const + { return m_indentationStr; } + private: // Number of column layers that a table layout has, default is 1; @@ -823,6 +882,8 @@ class PreparedTableLayout : public TableLayout // Numbers of lower most column that are visible size_t m_visibleLowermostColumnCount; + string m_indentationStr; + /** * @brief Recursive part of column layout preparation, see constructor documentation. * @param columns The list of columns to prepare. diff --git a/src/coreComponents/common/format/table/TableMpiComponents.cpp b/src/coreComponents/common/format/table/TableMpiComponents.cpp new file mode 100644 index 00000000000..ad553f72771 --- /dev/null +++ b/src/coreComponents/common/format/table/TableMpiComponents.cpp @@ -0,0 +1,160 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + + +/** + * @file TableMpiComponents.cpp + */ + +#include "TableMpiComponents.hpp" +#include "common/MpiWrapper.hpp" + +namespace geos +{ + +TableTextMpiOutput::TableTextMpiOutput( TableMpiLayout mpiLayout ): + TableTextFormatter(), + m_mpiLayout( mpiLayout ) +{} + +TableTextMpiOutput::TableTextMpiOutput( TableLayout const & tableLayout, + TableMpiLayout mpiLayout ): + TableTextFormatter( tableLayout ), + m_mpiLayout( mpiLayout ) +{} + +template<> +void TableTextMpiOutput::toStream< TableData >( std::ostream & tableOutput, + TableData const & tableData ) const +{ + TableTextMpiOutput::Status status { + // m_isMasterRank (only the master rank does the output of the header && bottom of the table) + MpiWrapper::commRank() == 0, + // m_isContributing (some ranks does not have any output to produce) + !tableData.getCellsData().empty(), + // m_hasContent + false, + // m_sepLine + "" + }; + + CellLayoutRows headerCellsLayout; + CellLayoutRows dataCellsLayout; + CellLayoutRows errorCellsLayout; + size_t tableTotalWidth = 0; + + { + ColumnWidthModifier const columnWidthModifier = [this, status]( stdVector< size_t > & columnsWidth ) { + stretchColumnsByRanks( columnsWidth, status ); + }; + initalizeTableGrids( m_tableLayout, tableData, + headerCellsLayout, dataCellsLayout, errorCellsLayout, + tableTotalWidth, columnWidthModifier ); + status.m_sepLine = string( tableTotalWidth, m_horizontalLine ); + } + + if( status.m_isMasterRank ) + { + outputTableHeader( tableOutput, m_tableLayout, headerCellsLayout, status.m_sepLine ); + tableOutput.flush(); + } + + outputTableDataToRank0( tableOutput, m_tableLayout, dataCellsLayout, status ); + + if( status.m_isMasterRank ) + { + outputTableFooter( tableOutput, m_tableLayout, errorCellsLayout, + status.m_sepLine, status.m_hasContent ); + tableOutput.flush(); + } +} + +void TableTextMpiOutput::stretchColumnsByRanks( stdVector< size_t > & columnsWidth, + TableTextMpiOutput::Status const & status ) const +{ + { // we ensure we have the correct amount of columns on all ranks (for correct MPI reduction operation) + size_t const rankColumnsCount = columnsWidth.size(); + size_t const maxRanksColumnsCount = MpiWrapper::max( rankColumnsCount ); + + // TODO: contribute to the new table error system with this one + if( status.m_isContributing ) + GEOS_ASSERT_EQ( rankColumnsCount, maxRanksColumnsCount ); + + columnsWidth.resize( maxRanksColumnsCount, 0 ); + } + + // the ranks that does not contribute must not interfere in the column width computing + if( !status.m_isContributing ) + std::fill( columnsWidth.begin(), columnsWidth.end(), 0 ); + + // we keep the largest column widths so we have the same layout on all ranks + MpiWrapper::allReduce( columnsWidth, columnsWidth, MpiWrapper::Reduction::Max ); +} + +void TableTextMpiOutput::outputTableDataToRank0( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & dataCellsLayout, + TableTextMpiOutput::Status & status ) const +{ + integer const ranksCount = MpiWrapper::commSize(); + + // master rank does the output directly to the output, other ranks will have to send it through a string. + std::ostringstream localStringStream; + std::ostream & rankOutput = status.m_isMasterRank ? tableOutput : localStringStream; + + if( status.m_isContributing ) + { + if( m_mpiLayout.m_separatorBetweenRanks ) + { + string const rankSepLine = GEOS_FMT( "{:-^{}}", m_mpiLayout.m_rankTitle, status.m_sepLine.size() - 2 ); + rankOutput << tableLayout.getIndentationStr() << m_verticalLine << rankSepLine << m_verticalLine << '\n'; + } + outputTableData( rankOutput, tableLayout, dataCellsLayout ); + } + + // all other ranks than rank 0 render their output in a string and comunicate its size + std::vector< integer > ranksStrsSizes = std::vector( ranksCount, 0 ); + string const rankStr = !status.m_isMasterRank && status.m_isContributing ? localStringStream.str() : ""; + integer const rankStrSize = rankStr.size(); + MpiWrapper::allgather( &rankStrSize, 1, ranksStrsSizes.data(), 1 ); + + // we compute the memory layout of the ranks strings + std::vector< integer > ranksStrsOffsets = std::vector( ranksCount, 0 ); + integer ranksStrsTotalSize = 0; + for( integer rankId = 1; rankId < ranksCount; ++rankId ) + { + ranksStrsOffsets[rankId] = ranksStrsTotalSize; + ranksStrsTotalSize += ranksStrsSizes[rankId]; + } + + // finally, we can send all text data to rank 0, then we output it in the output stream. + string ranksStrs = string( ranksStrsTotalSize, '\0' ); + MpiWrapper::gatherv( &rankStr[0], rankStrSize, + &ranksStrs[0], ranksStrsSizes.data(), ranksStrsOffsets.data(), + 0, MPI_COMM_GEOS ); + if( status.m_isMasterRank ) + { + for( integer rankId = 1; rankId < ranksCount; ++rankId ) + { + if( ranksStrsSizes[rankId] > 0 ) + { + status.m_hasContent = true; + tableOutput << string_view( &ranksStrs[ranksStrsOffsets[rankId]], ranksStrsSizes[rankId] ); + } + } + } +} + +} /* namespace geos */ diff --git a/src/coreComponents/common/format/table/TableMpiComponents.hpp b/src/coreComponents/common/format/table/TableMpiComponents.hpp new file mode 100644 index 00000000000..aeb37caa0ce --- /dev/null +++ b/src/coreComponents/common/format/table/TableMpiComponents.hpp @@ -0,0 +1,110 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file TableMpiComponents.hpp + * @brief contains variation of tables components for MPI communication. + */ + +#ifndef GEOS_COMMON_FORMAT_TABLE_TABLEMPICOMPONENTS_HPP +#define GEOS_COMMON_FORMAT_TABLE_TABLEMPICOMPONENTS_HPP + +#include "TableFormatter.hpp" + +namespace geos +{ + +/** + * @struct TableMpiLayout + * @brief Layout information specific to MPI distributed tables, completing those in TableLayout. + */ +struct TableMpiLayout +{ + /// Enable a separating line between ranks in the table output. + bool m_separatorBetweenRanks = false; + /// Title for each rank's section visible in the separating line. + string m_rankTitle; +}; + +/** + * @brief class to format data in a formatted text format, allowing contributions from multiple + * MPI ranks. + */ +class TableTextMpiOutput : public TableTextFormatter +{ +public: + /// base class + using Base = TableTextFormatter; + + /** + * @brief Construct a default Table Formatter without layout specification (to only insert data in it, + * without any column / title). Feature is not tested. + * @param mpiLayout MPI-specific layout information (default is having contiguous ranks data). + */ + TableTextMpiOutput( TableMpiLayout mpiLayout = TableMpiLayout() ); + + /** + * @brief Construct a new TableTextMpiOutput from a tableLayout + * @param tableLayout Contain all tableColumnData names and optionnaly the table title + * @param mpiLayout MPI-specific layout information (default is having contiguous ranks data). + */ + TableTextMpiOutput( TableLayout const & tableLayout, + TableMpiLayout mpiLayout = TableMpiLayout() ); + + /** + * @brief Convert a data source to a table string. + * @param tableData The data source to convert. + * @param outputStream The target output stream for rank 0, to output the table string representation + * of the TableData. Each rank contributing to the common rank 0 output stream + * with their local data. It may be the log or a file stream. + * @note This method must be called by all MPI ranks. + */ + template< typename DATASOURCE > + void toStream( std::ostream & outputStream, DATASOURCE const & tableData ) const; + +private: + + // hiding toString() methods as they are not implemented with MPI support. + using Base::toString; + + struct Status + { + bool const m_isMasterRank; + bool const m_isContributing; + bool m_hasContent; + string m_sepLine; + }; + + TableMpiLayout m_mpiLayout; + + /** + * @brief Expend the columns width to accomodate with the content of all MPI ranks. + * As it is based on MPI communications, every ranks must call this method. + * @param columnsWidth The array to store the resulting columns width in. + * @param tableGrid The grid of cells containing content. + */ + void stretchColumnsByRanks( stdVector< size_t > & columnsWidth, + Status const & status ) const; + + void outputTableDataToRank0( std::ostream & tableOutput, + PreparedTableLayout const & tableLayout, + CellLayoutRows const & dataCellsLayout, + Status & status ) const; + +}; + +} + +#endif /* GEOS_COMMON_FORMAT_TABLE_TABLEMPICOMPONENTS_HPP */ diff --git a/src/coreComponents/common/format/table/unitTests/CMakeLists.txt b/src/coreComponents/common/format/table/unitTests/CMakeLists.txt index f5a4a358e97..594f08577a7 100644 --- a/src/coreComponents/common/format/table/unitTests/CMakeLists.txt +++ b/src/coreComponents/common/format/table/unitTests/CMakeLists.txt @@ -2,6 +2,9 @@ set( gtest_geosx_tests testTable.cpp ) +set( gtest_geosx_mpi_tests + testMpiTable.cpp ) + set( dependencyList gtest common ${parallelDeps} ) # Add gtest C++ based tests @@ -16,3 +19,20 @@ foreach(test ${gtest_geosx_tests}) COMMAND ${test_name} ) endforeach() + +if( ENABLE_MPI ) + set( nranks 4 ) + + foreach( test ${gtest_geosx_mpi_tests} ) + get_filename_component( file_we ${test} NAME_WE ) + set( test_name ${file_we}_mpi ) + blt_add_executable( NAME ${test_name} + SOURCES ${test} + OUTPUT_DIR ${TEST_OUTPUT_DIRECTORY} + DEPENDS_ON ${dependencyList} ) + + geos_add_test( NAME ${test_name} + COMMAND ${test_name} -x ${nranks} + NUM_MPI_TASKS ${nranks} ) + endforeach() +endif() diff --git a/src/coreComponents/common/format/table/unitTests/testMpiTable.cpp b/src/coreComponents/common/format/table/unitTests/testMpiTable.cpp new file mode 100644 index 00000000000..7b01cf49515 --- /dev/null +++ b/src/coreComponents/common/format/table/unitTests/testMpiTable.cpp @@ -0,0 +1,149 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +// Source includes +#include "common/format/table/TableMpiComponents.hpp" +#include "common/initializeEnvironment.hpp" +#include "common/MpiWrapper.hpp" +// TPL includes +#include +#include + +using namespace geos; + +class MpiTestScope +{ +public: + + MpiTestScope( int argc, char * argv[] ) + { + ::testing::InitGoogleTest( &argc, argv ); + geos::setupEnvironment( argc, argv ); + } + + ~MpiTestScope() + { + geos::cleanupEnvironment(); + } + +}; + +TEST( testMpiTables, testDifferentRankData ) +{ + struct TestCase + { + stdVector< stdVector< std::pair< integer, real64 > > > m_ranksValues; + string m_expectedResult; + }; + + stdVector< TestCase > const testCases = + { + { + { // m_ranksValues: in this test, rank 2 has no value + { {1, 0.502} }, + { {2, 0.624}, {3, 0.791} }, + {}, + { {4, 0.243}, {5, 0.804}, {6, 0.302} }, + }, + "\n" // m_expectedResult + "-------------------------------------------\n" + "| Summary of negative pressure elements |\n" + "|-----------------------------------------|\n" + "| Global Id | pressure [Pa] |\n" + "|------------------|----------------------|\n" + "|------------Rank 0, 1 values-------------|\n" + "| 1 | 0.502 |\n" + "|------------Rank 1, 2 values-------------|\n" + "| 2 | 0.624 |\n" + "| 3 | 0.791 |\n" + "|------------Rank 3, 3 values-------------|\n" + "| 4 | 0.243 |\n" + "| 5 | 0.804 |\n" + "| 6 | 0.302 |\n" + "-------------------------------------------\n" + }, + { // m_ranksValues: in this test, rank 0 has no value + { + {}, + { {4, 0.243}, {5, 0.804}, {6, 0.302} }, + { {1, 0.502} }, + { {2, 0.624}, {3, 0.791} }, + }, + "\n" // m_expectedResult + "-------------------------------------------\n" + "| Summary of negative pressure elements |\n" + "|-----------------------------------------|\n" + "| Global Id | pressure [Pa] |\n" + "|------------------|----------------------|\n" + "|------------Rank 1, 3 values-------------|\n" + "| 4 | 0.243 |\n" + "| 5 | 0.804 |\n" + "| 6 | 0.302 |\n" + "|------------Rank 2, 1 values-------------|\n" + "| 1 | 0.502 |\n" + "|------------Rank 3, 2 values-------------|\n" + "| 2 | 0.624 |\n" + "| 3 | 0.791 |\n" + "-------------------------------------------\n" + }, + }; + for( TestCase const & testCase: testCases ) + { + int const rankId = MpiWrapper::commRank(); + int const nbRanks = MpiWrapper::commSize(); + if( nbRanks > 1 ) + { + ASSERT_EQ( nbRanks, 4 ); + + TableLayout const layout = TableLayout(). + setTitle( "Summary of negative pressure elements" ). + addColumns( { "Global Id", "pressure [Pa]" } ). + setDefaultHeaderAlignment( TableLayout::Alignment::left ); + TableData data; + auto const & rankTestData = testCase.m_ranksValues[rankId]; + + TableMpiLayout mpiLayout; + mpiLayout.m_separatorBetweenRanks = true; + + if( !rankTestData.empty() ) + { + mpiLayout.m_rankTitle = GEOS_FMT( "Rank {}, {} values", rankId, rankTestData.size() ); + for( auto const & [id, value] : rankTestData ) + { + data.addRow( id, value ); + } + } + + TableTextMpiOutput const formatter = TableTextMpiOutput( layout, mpiLayout ); + std::ostringstream oss; + formatter.toStream( oss, data ); + if( rankId == 0 ) + { + EXPECT_STREQ( testCase.m_expectedResult.data(), + oss.str().data() ); + } + } + } +} + +int main( int argc, char * * argv ) +{ + int r; + { + MpiTestScope testScope{ argc, argv }; + r = RUN_ALL_TESTS(); + } + return r; +} diff --git a/src/coreComponents/common/format/table/unitTests/testTable.cpp b/src/coreComponents/common/format/table/unitTests/testTable.cpp index a5201f10818..b8f76efa758 100644 --- a/src/coreComponents/common/format/table/unitTests/testTable.cpp +++ b/src/coreComponents/common/format/table/unitTests/testTable.cpp @@ -18,8 +18,6 @@ #include "common/format/table/TableFormatter.hpp" #include "common/format/table/TableLayout.hpp" #include "dataRepository/Group.hpp" -#include "common/DataTypes.hpp" -#include "LvArray/src/MallocBuffer.hpp" // TPL includes #include #include @@ -216,7 +214,7 @@ TEST( testTable, tableHiddenColumn ) TEST( testTable, tableMergeOverflowParadox ) { string const title = "Lorem Ipsum"; - TableLayout tableLayout( title, + TableLayout const tableLayout( title, { TableLayout::Column() .setName( "A" ), @@ -858,6 +856,28 @@ TEST( testTable, tableSpecialsValues ) ); } +TEST( testTable, testTitleWithNoColumnIndented ) +{ + TableLayout const tableLayout = TableLayout(). + setTitle( "Title" ). + setIndentation( 4 ). + setMargin( TableLayout::MarginValue::small ); + + TableData tableData; + tableData.addRow( "Global Id", 1234, 40, 5678, 60 ); + tableData.addRow( "pressure", 0.1234, 0.40, 0.5678, 0.60 ); + + TableTextFormatter const tableText( tableLayout ); + EXPECT_EQ( tableText.toString( tableData ), + "\n" + " -------------------------------------------\n" + " | Title |\n" + " |-----------------------------------------|\n" + " | Global Id | 1234 | 40 | 5678 | 60 |\n" + " | pressure | 0.1234 | 0.4 | 0.5678 | 0.6 |\n" + " -------------------------------------------\n" + ); +} int main( int argc, char * * argv ) { diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/parameters/KValueFlashParameters.cpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/parameters/KValueFlashParameters.cpp index 8a5e751d640..6e40b9848c5 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/parameters/KValueFlashParameters.cpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/parameters/KValueFlashParameters.cpp @@ -381,7 +381,7 @@ bool KValueFlashParameters< NUM_PHASE >::validateKValues( MultiFluidBase const * } hasAtLeastOneNegative = hasAtLeastOneNegative || hasNegative; hasAtLeastOneOneSided = hasAtLeastOneOneSided || (allMoreThanUnity || allLessThanUnity); - if( (allMoreThanUnity || allLessThanUnity || hasNegative) && tableData.getTableDataRows().size() < 5 ) + if( (allMoreThanUnity || allLessThanUnity || hasNegative) && tableData.getCellsData().size() < 5 ) { tableRow[0].value = phaseNames[phaseIndex+1]; tableRow[1].value = GEOS_FMT( "{0:.3e}", m_pressureValues[0][pressureIndex] ); @@ -397,7 +397,7 @@ bool KValueFlashParameters< NUM_PHASE >::validateKValues( MultiFluidBase const * } } - if( !tableData.getTableDataRows().empty()) + if( !tableData.getCellsData().empty()) { std::vector< TableLayout::Column > columns; columns.emplace_back( TableLayout::Column().setName( "Phase" ).setValuesAlignment( TableLayout::Alignment::left ) ); diff --git a/src/coreComponents/physicsSolvers/LogLevelsInfo.hpp b/src/coreComponents/physicsSolvers/LogLevelsInfo.hpp index e11aca92b2e..73037d8bb04 100644 --- a/src/coreComponents/physicsSolvers/LogLevelsInfo.hpp +++ b/src/coreComponents/physicsSolvers/LogLevelsInfo.hpp @@ -91,6 +91,12 @@ struct Solution static constexpr std::string_view getDescription() { return "Solution information (scaling, maximum changes, quality check)"; } }; +struct SolutionDetails +{ + static constexpr int getMinLogLevel() { return 2; } + static constexpr std::string_view getDescription() { return "Solution details (incoherent negative values ids)"; } +}; + struct SolverInitialization { static constexpr int getMinLogLevel() { return 1; } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CMakeLists.txt b/src/coreComponents/physicsSolvers/fluidFlow/CMakeLists.txt index 6eb55e9775a..1f6bd171acc 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CMakeLists.txt +++ b/src/coreComponents/physicsSolvers/fluidFlow/CMakeLists.txt @@ -39,9 +39,11 @@ set( fluidFlowSolvers_headers SinglePhaseHybridFVM.hpp SinglePhaseProppantBase.hpp StencilAccessors.hpp + SolutionCheckHelpers.hpp StencilDataCollection.hpp LogLevelsInfo.hpp kernels/MinPoreVolumeMaxPorosityKernel.hpp + kernels/SolutionCheckKernelsHelpers.hpp kernels/StencilWeightsUpdateKernel.hpp kernels/HybridFVMHelperKernels.hpp kernels/singlePhase/AccumulationKernels.hpp @@ -154,6 +156,7 @@ set( fluidFlowSolvers_sources SinglePhaseProppantBase.cpp SinglePhaseReactiveTransport.cpp SourceFluxStatistics.cpp + SolutionCheckHelpers.cpp StencilDataCollection.cpp kernels/singlePhase/proppant/ProppantFluxKernels.cpp kernels/compositional/AquiferBCKernel.cpp diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp index 40a8b4ddd62..89aec54eb4f 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp @@ -37,6 +37,7 @@ #include "physicsSolvers/LogLevelsInfo.hpp" #include "physicsSolvers/fluidFlow/FlowSolverBaseFields.hpp" #include "physicsSolvers/fluidFlow/CompositionalMultiphaseBaseFields.hpp" +#include "physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp" #include "physicsSolvers/fluidFlow/kernels/compositional/ResidualNormKernel.hpp" #include "physicsSolvers/fluidFlow/kernels/compositional/ThermalResidualNormKernel.hpp" #include "physicsSolvers/fluidFlow/kernels/compositional/SolutionScalingKernel.hpp" @@ -848,7 +849,11 @@ bool CompositionalMultiphaseFVM::checkSystemSolution( DomainPartition & domain, string const dofKey = dofManager.getKey( viewKeyStruct::elemDofFieldString() ); integer localCheck = 1; real64 minPres = 0.0, minDens = 0.0, minTotalDens = 0.0; - integer numNegPres = 0, numNegDens = 0, numNegTotalDens = 0; + integer numNegTotalDens = 0; + ElementsReporterBuffer rankNegPressureIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( getLogLevel() ) ? 16 : 0 }; + ElementsReporterBuffer rankNegDensityIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( this->getLogLevel() ) ? 16 : 0 }; forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &, MeshLevel & mesh, @@ -867,52 +872,57 @@ bool CompositionalMultiphaseFVM::checkSystemSolution( DomainPartition & domain, arrayView1d< real64 > pressureScalingFactor = subRegion.getField< flow::pressureScalingFactor >(); arrayView1d< real64 > temperatureScalingFactor = subRegion.getField< flow::temperatureScalingFactor >(); arrayView1d< real64 > compDensScalingFactor = subRegion.getField< flow::globalCompDensityScalingFactor >(); + auto const & cellLocalToGlobalIds = subRegion.localToGlobalMap(); + auto const negPresCollector = rankNegPressureIds.createCollector( cellLocalToGlobalIds ); + auto const negDensCollector = rankNegDensityIds.createCollector( cellLocalToGlobalIds ); + // check that pressure and component densities are non-negative // for thermal, check that temperature is above 273.15 K const integer temperatureOffset = m_numComponents+1; - auto const subRegionData = - m_isThermal - ? thermalCompositionalMultiphaseBaseKernels:: - SolutionCheckKernelFactory:: - createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, - m_allowNegativePressure, - m_scalingType, - scalingFactor, - pressure, - temperature, - compDens, - pressureScalingFactor, - temperatureScalingFactor, - compDensScalingFactor, - dofManager.rankOffset(), - m_numComponents, - dofKey, - subRegion, - localSolution, - temperatureOffset ) - : isothermalCompositionalMultiphaseBaseKernels:: - SolutionCheckKernelFactory:: - createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, - m_allowNegativePressure, - m_scalingType, - scalingFactor, - pressure, - compDens, - pressureScalingFactor, - compDensScalingFactor, - dofManager.rankOffset(), - m_numComponents, - dofKey, - subRegion, - localSolution ); + auto const subRegionData = m_isThermal ? + thermalCompositionalMultiphaseBaseKernels:: + SolutionCheckKernelFactory:: + createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, + m_allowNegativePressure, + m_scalingType, + scalingFactor, + pressure, + temperature, + compDens, + pressureScalingFactor, + temperatureScalingFactor, + compDensScalingFactor, + dofManager.rankOffset(), + m_numComponents, + dofKey, + subRegion, + localSolution, + negPresCollector, + negDensCollector, + temperatureOffset ) : + isothermalCompositionalMultiphaseBaseKernels:: + SolutionCheckKernelFactory:: + createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, + m_allowNegativePressure, + m_scalingType, + scalingFactor, + pressure, + compDens, + pressureScalingFactor, + compDensScalingFactor, + dofManager.rankOffset(), + m_numComponents, + dofKey, + subRegion, + localSolution, + negPresCollector, + negDensCollector ); localCheck = std::min( localCheck, subRegionData.localMinVal ); - minPres = std::min( minPres, subRegionData.localMinPres ); - minDens = std::min( minDens, subRegionData.localMinDens ); - minTotalDens = std::min( minTotalDens, subRegionData.localMinTotalDens ); - numNegPres += subRegionData.localNumNegPressures; - numNegDens += subRegionData.localNumNegDens; + minPres = std::min( minPres, subRegionData.localMinNegPres ); + minDens = std::min( minDens, subRegionData.localMinNegDens ); + minTotalDens = std::min( minTotalDens, subRegionData.localMinNegTotalDens ); numNegTotalDens += subRegionData.localNumNegTotalDens; } ); } ); @@ -920,23 +930,20 @@ bool CompositionalMultiphaseFVM::checkSystemSolution( DomainPartition & domain, minPres = MpiWrapper::min( minPres ); minDens = MpiWrapper::min( minDens ); minTotalDens = MpiWrapper::min( minTotalDens ); - numNegPres = MpiWrapper::sum( numNegPres ); - numNegDens = MpiWrapper::sum( numNegDens ); numNegTotalDens = MpiWrapper::sum( numNegTotalDens ); - if( numNegPres > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative pressure values: {}, minimum value: {} Pa", - getName(), numNegPres, fmt::format( "{:.{}f}", minPres, 3 ) ) ); - string const massUnit = m_useMass ? "kg/m3" : "mol/m3"; - if( numNegDens > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative component density values: {}, minimum value: {} {}}", - getName(), numNegDens, fmt::format( "{:.{}f}", minDens, 3 ), massUnit ) ); - if( minTotalDens > 0 ) + rankNegPressureIds.createOutput().outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative pressure", minPres, units::Unit::Pressure ); + + units::Unit const massUnit = m_useMass ? units::Unit::Density : units::Unit::MolarDensity; + rankNegDensityIds.createOutput().outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative component density", minDens, massUnit ); + if( numNegTotalDens > 0 ) + { GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative total density values: {}, minimum value: {} {}}", - getName(), minTotalDens, fmt::format( "{:.{}f}", minDens, 3 ), massUnit ) ); + GEOS_FMT( " {}: Number of negative total density values: {}, minimum value: {} {}", + getName(), numNegTotalDens, fmt::format( "{:.{}f}", minTotalDens, 3 ), units::getSymbol( massUnit ) ) ); + } return MpiWrapper::min( localCheck ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp index 68d14faaebd..ede2009eb62 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp @@ -622,7 +622,9 @@ bool CompositionalMultiphaseHybridFVM::checkSystemSolution( DomainPartition & do m_numComponents, elemDofKey, subRegion, - localSolution ); + localSolution, + ElementsReporterCollector::disabled(), + ElementsReporterCollector::disabled() ); localCheck = std::min( localCheck, subRegionData.localMinVal ); } ); diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp index 713367c05f0..5bd38614c34 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp @@ -39,6 +39,7 @@ #include "physicsSolvers/KernelLaunchSelectors.hpp" #include "physicsSolvers/fluidFlow/FlowSolverBaseFields.hpp" #include "physicsSolvers/fluidFlow/SinglePhaseBaseFields.hpp" +#include "physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp" #include "physicsSolvers/fluidFlow/kernels/singlePhase/MobilityKernel.hpp" #include "physicsSolvers/fluidFlow/kernels/singlePhase/SolutionCheckKernel.hpp" #include "physicsSolvers/fluidFlow/kernels/singlePhase/SolutionScalingKernel.hpp" @@ -1327,7 +1328,8 @@ bool SinglePhaseBase::checkSystemSolution( DomainPartition & domain, GEOS_MARK_FUNCTION; string const dofKey = dofManager.getKey( viewKeyStruct::elemDofFieldString() ); - integer numNegativePressures = 0; + ElementsReporterBuffer rankNegPressureIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( getLogLevel() ) ? 16 : 0 }; real64 minPressure = 0.0; forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &, @@ -1342,23 +1344,25 @@ bool SinglePhaseBase::checkSystemSolution( DomainPartition & domain, arrayView1d< globalIndex const > const dofNumber = subRegion.getReference< array1d< globalIndex > >( dofKey ); arrayView1d< integer const > const ghostRank = subRegion.ghostRank(); arrayView1d< real64 const > const pres = subRegion.getField< flow::pressure >(); - - auto const statistics = - singlePhaseBaseKernels::SolutionCheckKernel:: - launch< parallelDevicePolicy<> >( localSolution, rankOffset, dofNumber, ghostRank, pres, scalingFactor ); - - numNegativePressures += statistics.first; - minPressure = std::min( minPressure, statistics.second ); + auto const negPresCollector = rankNegPressureIds.createCollector( subRegion.localToGlobalMap().toViewConst() ); + + auto const statistics = singlePhaseBaseKernels::SolutionCheckKernel:: + launch< parallelDevicePolicy<> >( localSolution, + rankOffset, + dofNumber, + ghostRank, + pres, + scalingFactor, + negPresCollector ); + + minPressure = std::min( minPressure, statistics.minNegPres ); } ); } ); - numNegativePressures = MpiWrapper::sum( numNegativePressures ); - - if( numNegativePressures > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, GEOS_FMT( " {}: Number of negative pressure values: {}, minimum value: {} Pa", - getName(), numNegativePressures, fmt::format( "{:.{}f}", minPressure, 3 ) ) ); + rankNegPressureIds.createOutput().outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative pressure", minPressure, units::Unit::Pressure ); - return (m_allowNegativePressure || numNegativePressures == 0) ? 1 : 0; + return (m_allowNegativePressure || rankNegPressureIds.getSignaledElementsCount() == 0) ? 1 : 0; } void SinglePhaseBase::saveConvergedState( ElementSubRegionBase & subRegion ) const diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.cpp b/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.cpp new file mode 100644 index 00000000000..3d9366ef3b7 --- /dev/null +++ b/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.cpp @@ -0,0 +1,118 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file SolutionCheckKernel.cpp + */ + +#include "physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp" +#include "common/MpiWrapper.hpp" +#include "common/format/StringUtilities.hpp" +#include "common/format/table/TableMpiComponents.hpp" + +namespace geos +{ + +ElementsReporterBuffer::ElementsReporterBuffer( bool enabled, ElementCount maxCollectionSize ): + m_elementsCounter( enabled ? 1 : 0 ), + m_elementsBuffer( enabled ? maxCollectionSize : 0 ) +{ + if( enabled ) + { + m_elementsCounter.zero(); + m_elementsBuffer.zero(); + } +} + +ElementsReporterCollector +ElementsReporterBuffer::createCollector( arrayView1d< globalIndex const > const & localToGlobalId ) const +{ + return ElementsReporterCollector( m_elementsCounter, m_elementsBuffer, localToGlobalId ); +} + +ElementsReporterOutput ElementsReporterBuffer::createOutput() const +{ + m_elementsCounter.move( LvArray::MemorySpace::host, false ); + m_elementsBuffer.move( LvArray::MemorySpace::host, false ); + return ElementsReporterOutput( *this ); +} + + +ElementsReporterOutput::ElementsReporterOutput( ElementsReporterBuffer const & buffer ): + m_buffer( buffer ), + m_ranksSignaledElementsCount( MpiWrapper::sum( buffer.getSignaledElementsCount() ) ), + m_ranksCollectedElementsCount( MpiWrapper::sum( buffer.getCollectedElementsCount() ) ) +{} + +void ElementsReporterOutput::outputTooLowValues( string_view linesPrefix, + string_view valueNaming, + real64 minValue, + units::Unit unit ) const +{ + if( m_buffer.enabled() ) + { + if( m_ranksSignaledElementsCount > 0 ) + { + string const minValueStr = GEOS_FMT( "{:.{}f} [{}]", minValue, 3, units::getSymbol( unit ) ); + GEOS_LOG_RANK_0( GEOS_FMT( "{}{} {} values encountered. Minimum value: {}.", + linesPrefix, m_ranksSignaledElementsCount, valueNaming, minValueStr ) ); + + if( m_ranksCollectedElementsCount > 0 ) + { + TableLayout const layout = TableLayout(). + setTitle( GEOS_FMT( "Summary of {} elements", valueNaming ) ). + addColumns( { "Global Id", units::getDescription( unit ) } ). + enableLineBreak( false ). + setIndentation( linesPrefix.size() ). + setDefaultHeaderAlignment( TableLayout::Alignment::left ); + TableData data; + integer const signaledCount = m_buffer.getSignaledElementsCount(); + integer const collectedCount = m_buffer.getCollectedElementsCount(); + integer const omittedCount = signaledCount - collectedCount; + + TableMpiLayout mpiLayout; + mpiLayout.m_separatorBetweenRanks = true; + + if( signaledCount > 0 ) + { + mpiLayout.m_rankTitle = GEOS_FMT( "Rank {}, {} / {} values", + MpiWrapper::commRank(), collectedCount, signaledCount ); + + for( ElementReport const & report : m_buffer ) + { + data.addRow( report.m_id, report.m_value ); + } + + // adding one last line for signaling partial data & readability + if( omittedCount > 0 ) + { + data.addRow( "...", "..." ); + } + } + + TableTextMpiOutput const formatter = TableTextMpiOutput( layout, mpiLayout ); + formatter.toStream( std::cout, data ); + GEOS_LOG_RANK_0( '\n' ); + } + else + { + GEOS_LOG( GEOS_FMT( "{}Increase the log level to enable a reporting of the {} values.", + string( linesPrefix.size(), ' ' ), valueNaming ) ); + } + } + } +} + +} // namespace geos diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp b/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp new file mode 100644 index 00000000000..9965dcb3d75 --- /dev/null +++ b/src/coreComponents/physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp @@ -0,0 +1,190 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file SolutionCheckHelpers.hpp + */ + +#ifndef GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKHELPERS_HPP +#define GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKHELPERS_HPP + +#include "physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp" +#include "common/Units.hpp" + +namespace geos +{ + +/** + * @brief A class to report elements collected by the solver. + */ +class ElementsReporterOutput +{ +public: + + /// Type alias for elements count (e.g., localIndex, globalIndex). + using ElementCount = ElementsReporterCollector::ElementCount; + + /** + * @brief Construct a preallocated buffer for collecting element ids in kernels. + * @param buffer The buffer that will be utilized for counting & collecting elements IDs during kernel execution. + */ + ElementsReporterOutput( ElementsReporterBuffer const & buffer ); + + /** + * @return The number of ranks that have signaled an id. + */ + ElementCount getRanksSignaledIdsCount() const + { return m_ranksSignaledElementsCount; } + + /** + * @return The total count of collected elements across all ranks for signaling ids. + */ + ElementCount getRanksCollectedIdsCount() const + { return m_ranksCollectedElementsCount; } + + /** + * @brief Report elements with values below a specified threshold in the log: + * Outputs lines indicating which variables have collected element ids whose corresponding + * solution components are too low, potentially signaling underflow or numerical instability. + * @param linesPrefix Prefix for the line of text to be printed + * @param valueNaming The name used when referring to variables within this context (e.g., "pressure", "density"). + * @param minValue Minimum acceptable solution component values. Values below this threshold are reported. + * @param valueUnit Unit in which `minValue` is expressed. + */ + void outputTooLowValues( string_view linesPrefix, + string_view valueNaming, + real64 minValue, + units::Unit valueUnit ) const; + +private: + + /// Preallocated buffer for collecting ids. + ElementsReporterBuffer const & m_buffer; + + /// Count of signaled elements per rank. + ElementCount m_ranksSignaledElementsCount; + + /// Total collected signaling id count across ranks. + ElementCount m_ranksCollectedElementsCount; + +}; + +/** + * @brief A buffer to count and store element ids during kernel execution. + * This facilitates the reporting mechanism by allowing a preallocated space for storing & counting elements. + */ +class ElementsReporterBuffer +{ +public: + + /// Type alias for elements count (e.g., localIndex, globalIndex). + using ElementCount = ElementsReporterCollector::ElementCount; + + /** + * @brief Construct a preallocated buffer to collect a limited quantity of ids in kernels. + * @param maxCollectionSize Limit of the buffer. + * If 0, the buffering functionnality is disabled and only the counting is enabled. + */ + ElementsReporterBuffer( bool enabled, ElementCount maxCollectionSize ); + + /** + * @brief Transfers ownership of an ElementsReporterBuffer to another instance (move semantics). + */ + ElementsReporterBuffer( ElementsReporterBuffer && other ) = default; + + /** + * @brief Transfers ownership of an ElementsReporterBuffer to another instance (move semantics). + */ + ElementsReporterBuffer & operator=( ElementsReporterBuffer && other ) = default; + + /** + * @brief Copying prevented as it doesn't seem relevant / useful. + */ + ElementsReporterBuffer( ElementsReporterBuffer const & ) = delete; + + /** + * @brief Copying prevented as it doesn't seem relevant / useful. + */ + ElementsReporterBuffer & operator=( ElementsReporterBuffer const & ) = delete; + + /** + * @return the count of signaled elements. + */ + ElementCount getSignaledElementsCount() const + { return m_elementsCounter.empty() ? 0 : m_elementsCounter[0]; } + + /** + * @return the collected elements that could effectivly be stored (zero if no collection is enabled). + */ + ElementCount getCollectedElementsCount() const + { return LvArray::math::min( getSignaledElementsCount(), m_elementsBuffer.size() ); } + + /** + * @return a reference to an element report by its ID within the buffer (0 -> collected count-1). + */ + ElementReport const & operator[]( ElementCount id ) const + { return m_elementsBuffer[id]; } + + /** + * @return iterator pointing at beginning of collected elements in the buffer. + */ + auto begin() const + { return m_elementsBuffer.begin(); } + + /** + * @return iterator pointing after the last collected element in the buffer. + */ + auto end() const + { return m_elementsBuffer.begin() + getCollectedElementsCount(); } + + /** + * @return true when the collection of elements is enabled. + */ + bool enabled() const + { return !m_elementsCounter.empty(); } + + /** + * @return true when there are no elements collected (always false when enabled() is false). + */ + bool empty() const + { return getCollectedElementsCount() == 0; } + + /** + * @return true if the collection of elements completely fills the buffer. + */ + bool isComplete() const + { return getCollectedElementsCount() < getSignaledElementsCount(); } + + /** + * @return A view on the ids array owned by the instance. -> change comment to explain the interest for kernels + */ + ElementsReporterCollector createCollector( arrayView1d< globalIndex const > const & localToGlobalId ) const; + + ElementsReporterOutput createOutput() const; + +private: + + // array of one element to get benefit of managed host-device memory. + array1d< ElementCount > m_elementsCounter; + + // Preallocated array of ids of detected elements + array1d< ElementReport > m_elementsBuffer; + +}; + +} // namespace geos + + +#endif //GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKHELPERS_HPP diff --git a/src/coreComponents/physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp b/src/coreComponents/physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp new file mode 100644 index 00000000000..46a28a532b5 --- /dev/null +++ b/src/coreComponents/physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp @@ -0,0 +1,137 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 TotalEnergies + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file SolutionCheckKernelsHelpers.hpp + */ + +#ifndef GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKKERNELSHELPERS_HPP +#define GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKKERNELSHELPERS_HPP + +#include "common/DataTypes.hpp" + +namespace geos +{ + +class ElementsReporterBuffer; + +struct ElementReport +{ + /// the global id of the reported element + globalIndex m_id; + /// a value to report for the given element (i.e. a negative pressure, a density... Or if needed could be of a templated type for + /// composite values) + real64 m_value; +}; + +/** + * @brief Collects and reports elements ids and data using an atomic counter. + * This class provides functionality to collect data from multiple threads safely by incrementing + * through an atomic counter for each reported element's ID. The collected IDs are stored in a + * size limited buffer, which can be used later for reporting or analysis purposes. + */ +class ElementsReporterCollector +{ + friend class ElementsReporterBuffer; +public: + + using ElementCount = int32_t; + + /** + * @name Constructors + * @brief This object can be copied and moved as it only provides views to internal memory buffers. + */ + ///@{ + /** @cond DO_NOT_DOCUMENT */ + + ElementsReporterCollector( ElementsReporterCollector const & other ) = default; + + ElementsReporterCollector( ElementsReporterCollector && other ) = default; + + ElementsReporterCollector & operator=( ElementsReporterCollector const & other ) = default; + + ElementsReporterCollector & operator=( ElementsReporterCollector && other ) = default; + + /** @endcond */ + ///@} + + static ElementsReporterCollector disabled() + { + return ElementsReporterCollector( arrayView1d< ElementCount >(), + arrayView1d< ElementReport >(), + arrayView1d< globalIndex const >() ); + } + + /** + * @brief Collects a single element report and adds its ID to the output buffer if not disabled and + * there are available slots in the buffer. + * @tparam CollectorAtomicPolicy The atomic increment operation to use for thread-safe counter increments. + * @param report A constant reference to an `ElementReport` object containing data from a single element + */ + template< typename CollectorAtomicPolicy > + GEOS_HOST_DEVICE + void collectElement( CollectorAtomicPolicy, ElementReport const & report ) const + { + if( !m_elementsCounter.empty() ) + { + ElementCount const outputStart = RAJA::atomicInc< CollectorAtomicPolicy >( &m_elementsCounter[0] ); + + if( outputStart < m_elementsBuffer.size() ) + { + m_elementsBuffer[outputStart].m_id = m_localToGlobalId[report.m_id]; + m_elementsBuffer[outputStart].m_value = report.m_value; + } + } + } + + // // currently unused version for adding multiple ids from a given kernel + // template< typename AddedArray, typename ElementCount > + // GEOS_HOST_DEVICE + // void collectIds( AddedArray const & newIds, ElementCount newIdsCount ) + // { + // ElementCount const outputStart = RAJA::atomicAdd< CollectorAtomicPolicy >( &m_elementsCounter[0], newIdsCount ); + // ElementCount const maxNbIdToAdd = ElementCount( m_elementsBuffer.size() - outputStart ); + // newIdsCount = LvArray::math::min( newIdsCount, maxNbIdToAdd ); + // for( ElementCount i = 0; i < newIdsCount; ++i ) + // { + // m_elementsBuffer[outputStart + i] = newIds[i]; + // } + // } + +private: + + /// array of one element to get benefit of chai managed memory. + arrayView1d< ElementCount > m_elementsCounter; + + /// ids of detected elements, quantity limited to 'maxIdsCount' + arrayView1d< ElementReport > m_elementsBuffer; + + /// Maps local element IDs to their respective global indices. + arrayView1d< globalIndex const > m_localToGlobalId; + + ElementsReporterCollector( arrayView1d< ElementCount > const & elementsCounter, + arrayView1d< ElementReport > const & elementsBuffer, + arrayView1d< globalIndex const > const & localToGlobalId ): + m_elementsCounter( elementsCounter ), + m_elementsBuffer( elementsBuffer ), + m_localToGlobalId( localToGlobalId ) + {} + +}; + +} // namespace geos + + +#endif //GEOS_PHYSICSSOLVERS_FLUIDFLOW_SOLUTIONCHECKKERNELSHELPERS_HPP diff --git a/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/SolutionCheckKernel.hpp b/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/SolutionCheckKernel.hpp index 89d681bb02c..1e57ed629cf 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/SolutionCheckKernel.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/SolutionCheckKernel.hpp @@ -21,6 +21,7 @@ #define GEOS_PHYSICSSOLVERS_FLUIDFLOW_COMPOSITIONAL_SOLUTIONCHECKKERNEL_HPP #include "physicsSolvers/fluidFlow/kernels/compositional/SolutionScalingAndCheckingKernelBase.hpp" +#include "physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp" namespace geos { @@ -71,7 +72,9 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer integer const numComp, string const dofKey, ElementSubRegionBase const & subRegion, - arrayView1d< real64 const > const localSolution ) + arrayView1d< real64 const > const localSolution, + ElementsReporterCollector const & negPressureIds, + ElementsReporterCollector const & negDensityIds ) : Base( rankOffset, numComp, dofKey, @@ -84,44 +87,55 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer m_allowCompDensChopping( allowCompDensChopping ), m_allowNegativePressure( allowNegativePressure ), m_scalingFactor( scalingFactor ), - m_scalingType( scalingType ) + m_scalingType( scalingType ), + m_negPressureIds( negPressureIds ), + m_negDensityIds( negDensityIds ) {} /** - * @struct StackVariables + * @struct KernelStats * @brief Kernel variables located on the stack */ - struct StackVariables : public Base::StackVariables + struct KernelStats : public Base::StackVariables { GEOS_HOST_DEVICE - StackVariables() - { } - - StackVariables( real64 _localMinVal, - real64 _localMinPres, - real64 _localMinDens, - real64 _localMinTotalDens, - integer _localNumNegPressures, - integer _localNumNegDens, - integer _localNumNegTotalDens ) + KernelStats(): + Base::StackVariables() + {} + + KernelStats( real64 _localMinVal, + real64 _localNegMinPres, + real64 _localMinNegDens, + real64 _localMinNegTotalDens, + integer _localNumNegTotalDens ) : Base::StackVariables( _localMinVal ), - localMinPres( _localMinPres ), - localMinDens( _localMinDens ), - localMinTotalDens( _localMinTotalDens ), - localNumNegPressures( _localNumNegPressures ), - localNumNegDens( _localNumNegDens ), + localMinNegPres( _localNegMinPres ), + localMinNegDens( _localMinNegDens ), + localMinNegTotalDens( _localMinNegTotalDens ), localNumNegTotalDens( _localNumNegTotalDens ) { } - real64 localMinPres; - real64 localMinDens; - real64 localMinTotalDens; + real64 localMinNegPres; + real64 localMinNegDens; + real64 localMinNegTotalDens; - integer localNumNegPressures; - integer localNumNegDens; - integer localNumNegTotalDens; + localIndex localNumNegTotalDens; // Can only be 0 or 1 in each kernel + }; + + /** + * @struct StackVariables + * @brief Kernel variables located on the stack + */ + struct StackVariables : public KernelStats + { + GEOS_HOST_DEVICE + StackVariables(): + KernelStats() + { } + localIndex localNumNegPres; + localIndex localNumNegDens; }; /** @@ -132,19 +146,20 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer * @param[inout] kernelComponent the kernel component providing access to the compute function */ template< typename POLICY, typename KERNEL_TYPE > - static StackVariables + static KernelStats launch( localIndex const numElems, KERNEL_TYPE const & kernelComponent ) { - RAJA::ReduceMin< ReducePolicy< POLICY >, integer > globalMinVal( 1 ); + using reducePolicy = ReducePolicy< POLICY >; + using atomicPolicy = AtomicPolicy< POLICY >; - RAJA::ReduceMin< ReducePolicy< POLICY >, real64 > minPres( 0.0 ); - RAJA::ReduceMin< ReducePolicy< POLICY >, real64 > minDens( 0.0 ); - RAJA::ReduceMin< ReducePolicy< POLICY >, real64 > minTotalDens( 0.0 ); + RAJA::ReduceMin< reducePolicy, integer > globalMinVal( 1 ); - RAJA::ReduceSum< ReducePolicy< POLICY >, integer > numNegPressures( 0 ); - RAJA::ReduceSum< ReducePolicy< POLICY >, integer > numNegDens( 0 ); - RAJA::ReduceSum< ReducePolicy< POLICY >, integer > numNegTotalDens( 0 ); + RAJA::ReduceMin< reducePolicy, real64 > minPres( 0.0 ); + RAJA::ReduceMin< reducePolicy, real64 > minDens( 0.0 ); + RAJA::ReduceMin< reducePolicy, real64 > minTotalDens( 0.0 ); + + RAJA::ReduceSum< reducePolicy, localIndex > numNegTotalDens( 0 ); forAll< POLICY >( numElems, [=] GEOS_HOST_DEVICE ( localIndex const ei ) { @@ -153,28 +168,30 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer return; } - StackVariables stack; + StackVariables stack{}; kernelComponent.setup( ei, stack ); kernelComponent.compute( ei, stack ); globalMinVal.min( stack.localMinVal ); - minPres.min( stack.localMinPres ); - minDens.min( stack.localMinDens ); - minTotalDens.min( stack.localMinTotalDens ); + minPres.min( stack.localMinNegPres ); + minDens.min( stack.localMinNegDens ); + minTotalDens.min( stack.localMinNegTotalDens ); + + if( stack.localNumNegPres > 0 ) + kernelComponent.m_negPressureIds.collectElement( atomicPolicy{}, { ei, stack.localMinNegPres } ); + + if( stack.localNumNegDens > 0 ) + kernelComponent.m_negDensityIds.collectElement( atomicPolicy{}, { ei, stack.localMinNegDens } ); - numNegPressures += stack.localNumNegPressures; - numNegDens += stack.localNumNegDens; numNegTotalDens += stack.localNumNegTotalDens; } ); - return StackVariables( globalMinVal.get(), - minPres.get(), - minDens.get(), - minTotalDens.get(), - numNegPressures.get(), - numNegDens.get(), - numNegTotalDens.get() ); + return KernelStats( globalMinVal.get(), + minPres.get(), + minDens.get(), + minTotalDens.get(), + numNegTotalDens.get() ); } GEOS_HOST_DEVICE @@ -183,11 +200,11 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer { Base::setup( ei, stack ); - stack.localMinPres = 0.0; - stack.localMinDens = 0.0; - stack.localMinTotalDens = 0.0; + stack.localMinNegPres = 0.0; + stack.localMinNegDens = 0.0; + stack.localMinNegTotalDens = 0.0; - stack.localNumNegPressures = 0; + stack.localNumNegPres = 0; stack.localNumNegDens = 0; stack.localNumNegTotalDens = 0; } @@ -220,15 +237,15 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer bool const localScaling = m_scalingType == compositionalMultiphaseUtilities::ScalingType::Local; real64 const newPres = m_pressure[ei] + (localScaling ? m_pressureScalingFactor[ei] : m_scalingFactor) * m_localSolution[stack.localRow]; - if( newPres < 0 ) + if( newPres <= 0.0 ) { + stack.localNumNegPres = 1; + if( !m_allowNegativePressure ) - { stack.localMinVal = 0; - } - stack.localNumNegPressures += 1; - if( newPres < stack.localMinPres ) - stack.localMinPres = newPres; + + if( newPres < stack.localMinNegPres ) + stack.localMinNegPres = newPres; } // if component density chopping is not allowed, the time step fails if a component density is negative @@ -239,12 +256,13 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer for( integer ic = 0; ic < m_numComp; ++ic ) { real64 const newDens = m_compDens[ei][ic] + (localScaling ? m_compDensScalingFactor[ei] : m_scalingFactor) * m_localSolution[stack.localRow + ic + 1]; - if( newDens < 0 ) + if( newDens <= 0.0 ) { + stack.localNumNegDens = 1; stack.localMinVal = 0; - stack.localNumNegDens += 1; - if( newDens < stack.localMinDens ) - stack.localMinDens = newDens; + + if( newDens < stack.localMinNegDens ) + stack.localMinNegDens = newDens; } } } @@ -256,12 +274,13 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer real64 const newDens = m_compDens[ei][ic] + (localScaling ? m_compDensScalingFactor[ei] : m_scalingFactor) * m_localSolution[stack.localRow + ic + 1]; totalDens += ( newDens > 0.0 ) ? newDens : 0.0; } - if( totalDens < 0 ) + if( totalDens <= 0.0 ) { + stack.localNumNegTotalDens = 1; stack.localMinVal = 0; - stack.localNumNegTotalDens += 1; - if( totalDens < stack.localMinTotalDens ) - stack.localMinTotalDens = totalDens; + + if( totalDens < stack.localMinNegTotalDens ) + stack.localMinNegTotalDens = totalDens; } } @@ -282,6 +301,10 @@ class SolutionCheckKernel : public SolutionScalingAndCheckingKernelBase< integer /// scaling type (global or local) compositionalMultiphaseUtilities::ScalingType const m_scalingType; + ElementsReporterCollector const m_negPressureIds; + + ElementsReporterCollector const m_negDensityIds; + }; /** @@ -303,7 +326,7 @@ class SolutionCheckKernelFactory * @param[in] localSolution the Newton update */ template< typename POLICY > - static SolutionCheckKernel::StackVariables + static SolutionCheckKernel::KernelStats createAndLaunch( integer const allowCompDensChopping, integer const allowNegativePressure, compositionalMultiphaseUtilities::ScalingType const scalingType, @@ -316,11 +339,13 @@ class SolutionCheckKernelFactory integer const numComp, string const dofKey, ElementSubRegionBase & subRegion, - arrayView1d< real64 const > const localSolution ) + arrayView1d< real64 const > const localSolution, + ElementsReporterCollector const & negPressureIds, + ElementsReporterCollector const & negDensityIds ) { SolutionCheckKernel kernel( allowCompDensChopping, allowNegativePressure, scalingType, scalingFactor, pressure, compDens, pressureScalingFactor, compDensScalingFactor, rankOffset, - numComp, dofKey, subRegion, localSolution ); + numComp, dofKey, subRegion, localSolution, negPressureIds, negDensityIds ); return SolutionCheckKernel::launch< POLICY >( subRegion.size(), kernel ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/ThermalSolutionCheckKernel.hpp b/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/ThermalSolutionCheckKernel.hpp index b0447b33782..455eec757dd 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/ThermalSolutionCheckKernel.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/ThermalSolutionCheckKernel.hpp @@ -73,6 +73,8 @@ class SolutionCheckKernel : public isothermalCompositionalMultiphaseBaseKernels: string const dofKey, ElementSubRegionBase const & subRegion, arrayView1d< real64 const > const localSolution, + ElementsReporterCollector const & negPressureIds, + ElementsReporterCollector const & negDensityIds, integer const temperatureOffset ) : Base( allowCompDensChopping, allowNegativePressure, @@ -86,7 +88,9 @@ class SolutionCheckKernel : public isothermalCompositionalMultiphaseBaseKernels: numComp, dofKey, subRegion, - localSolution ), + localSolution, + negPressureIds, + negDensityIds ), m_temperature( temperature ), m_temperatureScalingFactor( temperatureScalingFactor ), m_temperatureOffset( temperatureOffset ) @@ -146,7 +150,7 @@ class SolutionCheckKernelFactory * @param[in] localSolution the Newton update */ template< typename POLICY > - static SolutionCheckKernel::StackVariables + static SolutionCheckKernel::KernelStats createAndLaunch( integer const allowCompDensChopping, integer const allowNegativePressure, compositionalMultiphaseUtilities::ScalingType const scalingType, @@ -162,11 +166,13 @@ class SolutionCheckKernelFactory string const dofKey, ElementSubRegionBase & subRegion, arrayView1d< real64 const > const localSolution, + ElementsReporterCollector const & negPressureIds, + ElementsReporterCollector const & negDensityIds, integer temperatureOffset ) { SolutionCheckKernel kernel( allowCompDensChopping, allowNegativePressure, scalingType, scalingFactor, pressure, temperature, compDens, pressureScalingFactor, compDensScalingFactor, temperatureScalingFactor, - rankOffset, numComp, dofKey, subRegion, localSolution, + rankOffset, numComp, dofKey, subRegion, localSolution, negPressureIds, negDensityIds, temperatureOffset ); return SolutionCheckKernel::launch< POLICY >( subRegion.size(), kernel ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/kernels/singlePhase/SolutionCheckKernel.hpp b/src/coreComponents/physicsSolvers/fluidFlow/kernels/singlePhase/SolutionCheckKernel.hpp index 24529dda836..16bd95ea27d 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/kernels/singlePhase/SolutionCheckKernel.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/kernels/singlePhase/SolutionCheckKernel.hpp @@ -22,6 +22,7 @@ #include "common/DataTypes.hpp" #include "common/GEOS_RAJA_Interface.hpp" +#include "physicsSolvers/fluidFlow/kernels/SolutionCheckKernelsHelpers.hpp" namespace geos { @@ -33,16 +34,24 @@ namespace singlePhaseBaseKernels struct SolutionCheckKernel { + + struct KernelStats + { + real64 minNegPres; + }; + template< typename POLICY > - static std::pair< integer, real64 > launch( arrayView1d< real64 const > const & localSolution, - globalIndex const rankOffset, - arrayView1d< globalIndex const > const & dofNumber, - arrayView1d< integer const > const & ghostRank, - arrayView1d< real64 const > const & pres, - real64 const scalingFactor ) + static KernelStats launch( arrayView1d< real64 const > const & localSolution, + globalIndex const rankOffset, + arrayView1d< globalIndex const > const & dofNumber, + arrayView1d< integer const > const & ghostRank, + arrayView1d< real64 const > const & pres, + real64 const scalingFactor, + ElementsReporterCollector const & negPressureIds ) { - RAJA::ReduceSum< ReducePolicy< POLICY >, integer > numNegativePressures( 0 ); - RAJA::ReduceMin< ReducePolicy< POLICY >, real64 > minPres( 0.0 ); + using reducePolicy = ReducePolicy< POLICY >; + using atomicPolicy = AtomicPolicy< POLICY >; + RAJA::ReduceMin< reducePolicy, real64 > minNegPres( 0.0 ); forAll< POLICY >( dofNumber.size(), [=] GEOS_HOST_DEVICE ( localIndex const ei ) { @@ -51,16 +60,16 @@ struct SolutionCheckKernel localIndex const lid = dofNumber[ei] - rankOffset; real64 const newPres = pres[ei] + scalingFactor * localSolution[lid]; - if( newPres < 0.0 ) + if( newPres <= 0.0 ) { - numNegativePressures += 1; - minPres.min( newPres ); + minNegPres.min( newPres ); + negPressureIds.collectElement( atomicPolicy{}, { ei, newPres } ); } } } ); - return { numNegativePressures.get(), minPres.get() }; + return KernelStats{ minNegPres.get() }; } }; diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp index 9f18e9741da..a786a4841e8 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/CompositionalMultiphaseWell.cpp @@ -32,6 +32,7 @@ #include "mesh/WellElementSubRegion.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" #include "physicsSolvers/LogLevelsInfo.hpp" +#include "physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp" #include "physicsSolvers/fluidFlow/wells/LogLevelsInfo.hpp" #include "physicsSolvers/fluidFlow/wells/WellSolverBaseFields.hpp" #include "physicsSolvers/fluidFlow/wells/WellFields.hpp" @@ -1614,7 +1615,11 @@ CompositionalMultiphaseWell::checkSystemSolution( DomainPartition & domain, string const wellDofKey = dofManager.getKey( wellElementDofName() ); integer localCheck = 1; real64 minPres = 0.0, minDens = 0.0, minTotalDens = 0.0; - integer numNegPres = 0, numNegDens = 0, numNegTotalDens = 0; + integer numNegTotalDens = 0; + ElementsReporterBuffer rankNegPressureIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( getLogLevel() ) ? 16 : 0 }; + ElementsReporterBuffer rankNegDensityIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( this->getLogLevel() ) ? 16 : 0 }; forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &, MeshLevel & mesh, @@ -1636,53 +1641,57 @@ CompositionalMultiphaseWell::checkSystemSolution( DomainPartition & domain, arrayView1d< real64 > pressureScalingFactor = subRegion.getField< well::pressureScalingFactor >(); arrayView1d< real64 > temperatureScalingFactor = subRegion.getField< well::temperatureScalingFactor >(); arrayView1d< real64 > compDensScalingFactor = subRegion.getField< well::globalCompDensityScalingFactor >(); + auto const & cellLocalToGlobalIds = subRegion.localToGlobalMap(); + auto const negPresCollector = rankNegPressureIds.createCollector( cellLocalToGlobalIds ); + auto const negDensCollector = rankNegDensityIds.createCollector( cellLocalToGlobalIds ); // check that pressure and component densities are non-negative // for thermal, check that temperature is above 273.15 K const integer temperatureOffset = m_numComponents+2; - auto const subRegionData = - m_isThermal - ? thermalCompositionalMultiphaseBaseKernels:: - SolutionCheckKernelFactory:: - createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, - m_allowNegativePressure, - m_scalingType, - scalingFactor, - pressure, - temperature, - compDens, - pressureScalingFactor, - temperatureScalingFactor, - compDensScalingFactor, - dofManager.rankOffset(), - m_numComponents, - wellDofKey, - subRegion, - localSolution, - temperatureOffset ) - : isothermalCompositionalMultiphaseBaseKernels:: - SolutionCheckKernelFactory:: - createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, - m_allowNegativePressure, - m_scalingType, - scalingFactor, - pressure, - compDens, - pressureScalingFactor, - compDensScalingFactor, - dofManager.rankOffset(), - m_numComponents, - wellDofKey, - subRegion, - localSolution ); + auto const subRegionData = m_isThermal ? + thermalCompositionalMultiphaseBaseKernels:: + SolutionCheckKernelFactory:: + createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, + m_allowNegativePressure, + m_scalingType, + scalingFactor, + pressure, + temperature, + compDens, + pressureScalingFactor, + temperatureScalingFactor, + compDensScalingFactor, + dofManager.rankOffset(), + m_numComponents, + wellDofKey, + subRegion, + localSolution, + negPresCollector, + negDensCollector, + temperatureOffset ) : + isothermalCompositionalMultiphaseBaseKernels:: + SolutionCheckKernelFactory:: + createAndLaunch< parallelDevicePolicy<> >( m_allowCompDensChopping, + m_allowNegativePressure, + m_scalingType, + scalingFactor, + pressure, + compDens, + pressureScalingFactor, + compDensScalingFactor, + dofManager.rankOffset(), + m_numComponents, + wellDofKey, + subRegion, + localSolution, + negPresCollector, + negDensCollector ); localCheck = std::min( localCheck, subRegionData.localMinVal ); - minPres = std::min( minPres, subRegionData.localMinPres ); - minDens = std::min( minDens, subRegionData.localMinDens ); - minTotalDens = std::min( minTotalDens, subRegionData.localMinTotalDens ); - numNegPres += subRegionData.localNumNegPressures; - numNegDens += subRegionData.localNumNegDens; + minPres = std::min( minPres, subRegionData.localMinNegPres ); + minDens = std::min( minDens, subRegionData.localMinNegDens ); + minTotalDens = std::min( minTotalDens, subRegionData.localMinNegTotalDens ); numNegTotalDens += subRegionData.localNumNegTotalDens; } ); } ); @@ -1690,23 +1699,20 @@ CompositionalMultiphaseWell::checkSystemSolution( DomainPartition & domain, minPres = MpiWrapper::min( minPres ); minDens = MpiWrapper::min( minDens ); minTotalDens = MpiWrapper::min( minTotalDens ); - numNegPres = MpiWrapper::sum( numNegPres ); - numNegDens = MpiWrapper::sum( numNegDens ); numNegTotalDens = MpiWrapper::sum( numNegTotalDens ); - if( numNegPres > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative well pressure values: {}, minimum value: {} Pa", - getName(), numNegPres, fmt::format( "{:.{}f}", minPres, 3 ) ) ); - string const massUnit = m_useMass ? "kg/m3" : "mol/m3"; - if( numNegDens > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative well component density values: {}, minimum value: {} {} ", - getName(), numNegDens, fmt::format( "{:.{}f}", minDens, 3 ), massUnit ) ); - if( minTotalDens > 0 ) - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative total well density values: {}, minimum value: {} {} ", - getName(), minTotalDens, fmt::format( "{:.{}f}", minDens, 3 ), massUnit ) ); + rankNegPressureIds.createOutput().outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative pressure", minPres, units::Unit::Pressure ); + + units::Unit const massUnit = m_useMass ? units::Unit::Density : units::Unit::MolarDensity; + rankNegDensityIds.createOutput().outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative component density", minDens, massUnit ); + if( numNegTotalDens > 0 ) + { + GEOS_LOG_LEVEL_RANK_0( logInfo::SolutionDetails, + GEOS_FMT( " {}: Number of negative total density values: {}, minimum value: {} {}", + getName(), numNegTotalDens, fmt::format( "{:.{}f}", minTotalDens, 3 ), units::getSymbol( massUnit ) ) ); + } return MpiWrapper::min( localCheck ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp index 6414693bf22..875ddfba049 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/wells/SinglePhaseWell.cpp @@ -33,6 +33,7 @@ #include "physicsSolvers/LogLevelsInfo.hpp" #include "physicsSolvers/fluidFlow/wells/LogLevelsInfo.hpp" #include "physicsSolvers/fluidFlow/SinglePhaseBase.hpp" +#include "physicsSolvers/fluidFlow/SolutionCheckHelpers.hpp" #include "physicsSolvers/fluidFlow/FlowSolverBaseFields.hpp" #include "physicsSolvers/fluidFlow/wells/WellFields.hpp" #include "physicsSolvers/fluidFlow/wells/WellSolverBaseFields.hpp" @@ -1016,8 +1017,9 @@ bool SinglePhaseWell::checkSystemSolution( DomainPartition & domain, GEOS_MARK_FUNCTION; string const wellDofKey = dofManager.getKey( wellElementDofName() ); - integer numNegativePressures = 0; - real64 minPressure = 0.0; + ElementsReporterBuffer rankNegPressureIds{ isLogLevelActive< logInfo::Solution >( getLogLevel() ), + isLogLevelActive< logInfo::SolutionDetails >( getLogLevel() ) ? 16 : 0 }; + real64 minNegPres = 0.0; forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&] ( string const &, MeshLevel const & mesh, @@ -1041,25 +1043,26 @@ bool SinglePhaseWell::checkSystemSolution( DomainPartition & domain, arrayView1d< real64 const > const & pres = subRegion.getField< well::pressure >(); - auto const statistics = - singlePhaseBaseKernels::SolutionCheckKernel:: - launch< parallelDevicePolicy<> >( localSolution, rankOffset, dofNumber, ghostRank, pres, scalingFactor ); + auto const negPresCollector = rankNegPressureIds.createCollector( subRegion.localToGlobalMap().toViewConst() ); - numNegativePressures += statistics.first; - minPressure = std::min( minPressure, statistics.second ); + auto const results = singlePhaseBaseKernels::SolutionCheckKernel:: + launch< parallelDevicePolicy<> >( localSolution, + rankOffset, + dofNumber, + ghostRank, + pres, + scalingFactor, + negPresCollector ); + + minNegPres = std::min( minNegPres, results.minNegPres ); } ); } ); - numNegativePressures = MpiWrapper::sum( numNegativePressures ); - - if( numNegativePressures > 0 ) - { - GEOS_LOG_LEVEL_RANK_0( logInfo::Solution, - GEOS_FMT( " {}: Number of negative pressure values: {}, minimum value: {} Pa", - getName(), numNegativePressures, fmt::format( "{:.{}f}", minPressure, 3 ) ) ); - } + ElementsReporterOutput const rankNegPressureIdsOutput = rankNegPressureIds.createOutput(); + rankNegPressureIdsOutput.outputTooLowValues( GEOS_FMT( " {}: ", getName() ), + "negative pressure", minNegPres, units::Unit::Pressure ); - return (m_allowNegativePressure || numNegativePressures == 0) ? 1 : 0; + return (m_allowNegativePressure || rankNegPressureIdsOutput.getRanksSignaledIdsCount() == 0) ? 1 : 0; } void