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/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/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/mesh/ElementRegionManager.cpp b/src/coreComponents/mesh/ElementRegionManager.cpp index 97b5eeca57f..01d4f44ae88 100644 --- a/src/coreComponents/mesh/ElementRegionManager.cpp +++ b/src/coreComponents/mesh/ElementRegionManager.cpp @@ -14,12 +14,22 @@ */ #include +#include #include #include "ElementRegionManager.hpp" +#include "codingUtilities/RTTypes.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" #include "mesh/mpiCommunications/CommunicationTools.hpp" #include "SurfaceElementRegion.hpp" #include "constitutive/ConstitutiveManager.hpp" @@ -29,6 +39,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 @@ -228,6 +240,82 @@ void ElementRegionManager::generateWells( CellBlockManagerABC const & cellBlockM } ); + int rank = MpiWrapper::commRank(); + // int size = MpiWrapper::commSize(); + + 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 ) + { + 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 > perfLocation = wellSubRegion.getPerforationData()->getLocation(); + + 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 ); + + 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 ); + } + } + + + if( rank == 0 ) + { + std::cout << "destBuffer "<< rcvPerfo[0]<< " "<getWellElements()[iperf], + perfLocation[iperf], + resultRegionName[0], + setSubRegionName[0], + cellId ); + } + + } + } + 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 nodeManager.setMaxGlobalIndex(); } diff --git a/src/coreComponents/mesh/generators/WellGeneratorBase.cpp b/src/coreComponents/mesh/generators/WellGeneratorBase.cpp index 0fe4dc833a0..40da7884cb4 100644 --- a/src/coreComponents/mesh/generators/WellGeneratorBase.cpp +++ b/src/coreComponents/mesh/generators/WellGeneratorBase.cpp @@ -140,7 +140,6 @@ void WellGeneratorBase::generateWellGeometry( ) if( isLogLevelActive< logInfo::GenerateWell >( this->getLogLevel() ) && MpiWrapper::commRank() == 0 ) { logInternalWell(); - logPerforationTable(); } } @@ -557,18 +556,4 @@ void WellGeneratorBase::logInternalWell() const GEOS_LOG_RANK_0( wellFormatter.toString( wellData )); } -void WellGeneratorBase::logPerforationTable() const -{ - TableData dataPerforation; - for( globalIndex iperf = 0; iperf < m_numPerforations; ++iperf ) - { - dataPerforation.addRow( iperf, m_perfCoords[iperf], m_perfElemId[iperf] ); - } - - TableLayout const layoutPerforation ( GEOS_FMT( "Well '{}' Perforation Table", getName()), - { "Perforation no.", "Coordinates", "Well element no." } ); - TableTextFormatter const logPerforation( layoutPerforation ); - GEOS_LOG_RANK_0( logPerforation.toString( dataPerforation )); -} - } diff --git a/src/coreComponents/mesh/generators/WellGeneratorBase.hpp b/src/coreComponents/mesh/generators/WellGeneratorBase.hpp index e9a2dd81046..0367b2f0109 100644 --- a/src/coreComponents/mesh/generators/WellGeneratorBase.hpp +++ b/src/coreComponents/mesh/generators/WellGeneratorBase.hpp @@ -306,7 +306,6 @@ class WellGeneratorBase : public MeshComponentBase /// @cond DO_NOT_DOCUMENT void logInternalWell() const; - void logPerforationTable() const; /// @endcond /// Global number of perforations 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