diff --git a/src/coreComponents/codingUtilities/UnitTestUtilities.hpp b/src/coreComponents/codingUtilities/UnitTestUtilities.hpp index 3755134d630..6eddbf6e1ce 100644 --- a/src/coreComponents/codingUtilities/UnitTestUtilities.hpp +++ b/src/coreComponents/codingUtilities/UnitTestUtilities.hpp @@ -82,8 +82,8 @@ T expected( T expectedSerial, constexpr real64 DEFAULT_ABS_TOL = 1E-12; constexpr real64 DEFAULT_REL_TOL = std::numeric_limits< real64 >::epsilon(); -::testing::AssertionResult checkRelativeErrorFormat( const char *, const char *, const char *, const char *, - real64 const v1, real64 const v2, real64 const relTol, real64 const absTol ) +inline ::testing::AssertionResult checkRelativeErrorFormat( const char *, const char *, const char *, const char *, + real64 const v1, real64 const v2, real64 const relTol, real64 const absTol ) { real64 const delta = std::abs( v1 - v2 ); real64 const value = std::max( std::abs( v1 ), std::abs( v2 ) ); @@ -99,23 +99,23 @@ ::testing::AssertionResult checkRelativeErrorFormat( const char *, const char *, return ::testing::AssertionSuccess(); } -::testing::AssertionResult checkRelativeErrorFormat( char const * a, char const * b, char const * c, - real64 const v1, real64 const v2, real64 const relTol ) +inline ::testing::AssertionResult checkRelativeErrorFormat( char const * a, char const * b, char const * c, + real64 const v1, real64 const v2, real64 const relTol ) { return checkRelativeErrorFormat( a, b, c, "DEFAULT_ABS_TOL", v1, v2, relTol, DEFAULT_ABS_TOL ); } -void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol ) +inline void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol ) { EXPECT_PRED_FORMAT3( checkRelativeErrorFormat, v1, v2, relTol ); } -void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, real64 const absTol ) +inline void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, real64 const absTol ) { EXPECT_PRED_FORMAT4( checkRelativeErrorFormat, v1, v2, relTol, absTol ); } -void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, string const & name ) +inline void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, string const & name ) { SCOPED_TRACE( name ); EXPECT_PRED_FORMAT3( checkRelativeErrorFormat, v1, v2, relTol ); } -void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, real64 const absTol, string const & name ) +inline void checkRelativeError( real64 const v1, real64 const v2, real64 const relTol, real64 const absTol, string const & name ) { SCOPED_TRACE( name ); EXPECT_PRED_FORMAT4( checkRelativeErrorFormat, v1, v2, relTol, absTol ); diff --git a/src/coreComponents/common/DataTypes.cpp b/src/coreComponents/common/DataTypes.cpp index f2de87138e6..c0231da13d8 100644 --- a/src/coreComponents/common/DataTypes.cpp +++ b/src/coreComponents/common/DataTypes.cpp @@ -185,8 +185,8 @@ rtTypes::RegexMapType rtTypes::createBasicTypesRegexMap() string_view const groupNameDesc = "Input value must be a string that cannot be empty and contains only upper/lower letters, digits, and the characters . - _"; string_view const groupNameRegex = "[a-zA-Z0-9.\\-_]+"; - string_view const groupNameRefDesc = "Input value must be a string that can contain only upper/lower letters, digits, and the characters . - _ /"; - string_view const groupNameRefRegex = "[a-zA-Z0-9.\\-_/]*"; + string_view const groupNameRefDesc = "Input value must be a string that can contain only upper/lower letters, digits, and the characters . - _ / *"; + string_view const groupNameRefRegex = "[a-zA-Z0-9.\\-_/*]*"; // Build master list of regexes diff --git a/src/coreComponents/common/Units.hpp b/src/coreComponents/common/Units.hpp index d760217fbb8..e072dfe460b 100644 --- a/src/coreComponents/common/Units.hpp +++ b/src/coreComponents/common/Units.hpp @@ -82,6 +82,18 @@ enum Unit : integer /// Solubility in g/L Solubility, + + /// Mass in kg + Mass, + + /// Mole in mol + Mole, + + /// Mass rate in kg/s + MassRate, + + /// Mole rate in mol/s + MoleRate, }; @@ -104,6 +116,10 @@ constexpr inline std::string_view getDescription( Unit unit ) case Enthalpy: return "enthalpy [J/kg]"; case Density: return "density [kg/m3]"; case Solubility: return "solubility [g/L]"; + case Mass: return "mass [kg]"; + case Mole: return "mole [mol]"; + case MassRate: return "mass rate [kg/s]"; + case MoleRate: return "mole rate [mol/s]"; } } @@ -122,6 +138,14 @@ constexpr inline std::string_view getSymbol( Unit unit ) case TemperatureInC: return "C"; case Distance: return "m"; case Time: return "s"; + case Viscosity: return "Pa*s"; + case Enthalpy: return "J/kg"; + case Density: return "kg/m3"; + case Solubility: return "g/L"; + case Mass: return "kg"; + case Mole: return "mol"; + case MassRate: return "kg/s"; + case MoleRate: return "mol/s"; } } @@ -143,6 +167,14 @@ inline string formatValue( real64 value, Unit unit ) case TemperatureInC: return GEOS_FMT( "temperature of {} [K]", convertCToK( value ) ); case Distance: return GEOS_FMT( "distance of {} [s]", value ); case Time: return GEOS_FMT( "time of {} [s]", value ); + case Viscosity: return GEOS_FMT( "viscosity of {} [Pa*s]", value ); + case Enthalpy: return GEOS_FMT( "enthalpy of {} [J/kg]", value ); + case Density: return GEOS_FMT( "density of {} [kg/m3]", value ); + case Solubility: return GEOS_FMT( "solubility of {} [g/L]", value ); + case Mass: return GEOS_FMT( "mass of {} [kg]", value ); + case Mole: return GEOS_FMT( "mole of {} [mol]", value ); + case MassRate: return GEOS_FMT( "mass rate of {} [kg/s]", value ); + case MoleRate: return GEOS_FMT( "mole rate of {} [mol/s]", value ); } } diff --git a/src/coreComponents/dataRepository/Group.cpp b/src/coreComponents/dataRepository/Group.cpp index 5be399e8449..3958a2313c8 100644 --- a/src/coreComponents/dataRepository/Group.cpp +++ b/src/coreComponents/dataRepository/Group.cpp @@ -644,9 +644,8 @@ void Group::postRestartInitializationRecursive() void Group::enableLogLevelInput() { - string const logLevelString = "logLevel"; - - registerWrapper( logLevelString, &m_logLevel ). + // TODO : Improve the Log Level description to clearly assign a usecase per log level (incoming PR). + registerWrapper( viewKeyStruct::logLevelString(), &m_logLevel ). setApplyDefaultValue( 0 ). setInputFlag( InputFlags::OPTIONAL ). setDescription( "Log level" ); diff --git a/src/coreComponents/dataRepository/Group.hpp b/src/coreComponents/dataRepository/Group.hpp index 6192f3cddab..84ad45f1fbe 100644 --- a/src/coreComponents/dataRepository/Group.hpp +++ b/src/coreComponents/dataRepository/Group.hpp @@ -1414,6 +1414,15 @@ class Group */ void setInputFlags( InputFlags flags ) { m_input_flags = flags; } + /** + * @brief Structure to hold scoped key names + */ + struct viewKeyStruct + { + /// @return String for the logLevel wrapper + static constexpr char const * logLevelString() { return "logLevel"; } + }; + ///@} /** diff --git a/src/coreComponents/dataRepository/Wrapper.hpp b/src/coreComponents/dataRepository/Wrapper.hpp index 404225490af..b3193af0ee5 100644 --- a/src/coreComponents/dataRepository/Wrapper.hpp +++ b/src/coreComponents/dataRepository/Wrapper.hpp @@ -893,6 +893,15 @@ class Wrapper final : public WrapperBase return *this; } + /** + * @copydoc WrapperBase::appendDescription(string const &) + */ + Wrapper< T > & appendDescription( string const & description ) + { + WrapperBase::appendDescription( description ); + return *this; + } + /** * @copydoc WrapperBase::setRegisteringObjects(string const &) */ diff --git a/src/coreComponents/dataRepository/WrapperBase.hpp b/src/coreComponents/dataRepository/WrapperBase.hpp index 542863e03b9..d44e24c2e6b 100644 --- a/src/coreComponents/dataRepository/WrapperBase.hpp +++ b/src/coreComponents/dataRepository/WrapperBase.hpp @@ -508,6 +508,17 @@ class WrapperBase return *this; } + /** + * @brief Add up more text to the existing description string of the wrapper. + * @param description the description to add to the end of the previous one. + * @return a pointer to this wrapper + */ + WrapperBase & appendDescription( string const & description ) + { + m_description += description; + return *this; + } + /** * @brief Get the description string of the wrapper. * @return this wrapper's description string diff --git a/src/coreComponents/fieldSpecification/CMakeLists.txt b/src/coreComponents/fieldSpecification/CMakeLists.txt index 10ad2e8ec1a..ae31a4157fc 100644 --- a/src/coreComponents/fieldSpecification/CMakeLists.txt +++ b/src/coreComponents/fieldSpecification/CMakeLists.txt @@ -10,6 +10,7 @@ set( fieldSpecification_headers TractionBoundaryCondition.hpp AquiferBoundaryCondition.hpp PerfectlyMatchedLayer.hpp + SourceFluxStatistics.hpp ) # @@ -24,6 +25,7 @@ set( fieldSpecification_sources TractionBoundaryCondition.cpp AquiferBoundaryCondition.cpp PerfectlyMatchedLayer.cpp + SourceFluxStatistics.cpp ) set( dependencyList ${parallelDeps} functions linearAlgebra ) diff --git a/src/coreComponents/fieldSpecification/SourceFluxBoundaryCondition.cpp b/src/coreComponents/fieldSpecification/SourceFluxBoundaryCondition.cpp index acc3cc39981..53007480a30 100644 --- a/src/coreComponents/fieldSpecification/SourceFluxBoundaryCondition.cpp +++ b/src/coreComponents/fieldSpecification/SourceFluxBoundaryCondition.cpp @@ -18,6 +18,7 @@ */ #include "SourceFluxBoundaryCondition.hpp" +#include "physicsSolvers/fluidFlow/CompositionalMultiphaseBase.hpp" namespace geos { @@ -29,6 +30,18 @@ SourceFluxBoundaryCondition::SourceFluxBoundaryCondition( string const & name, G getWrapper< string >( FieldSpecificationBase::viewKeyStruct::fieldNameString() ). setInputFlag( InputFlags::FALSE ); setFieldName( catalogName() ); + + getWrapper< string >( FieldSpecificationBase::viewKeyStruct::functionNameString() ). + setDescription( GEOS_FMT( "Name of a function that specifies the variation of the production rate variations of this {}." + "Multiplied by {}. If no function is provided, a constant value of 1 is used." + "The producted fluid rate unit is in kg by default, or in mole if the flow solver has a {} of 0.", + catalogName(), + FieldSpecificationBase::viewKeyStruct::scaleString(), + CompositionalMultiphaseBase::viewKeyStruct::useMassFlagString() ) ); + + getWrapper< real64 >( FieldSpecificationBase::viewKeyStruct::scaleString() ). + setDescription( GEOS_FMT( "Multiplier of the {0} value. If no {0} is provided, this value is used directly.", + FieldSpecificationBase::viewKeyStruct::functionNameString() ) ); } REGISTER_CATALOG_ENTRY( FieldSpecificationBase, SourceFluxBoundaryCondition, string const &, Group * const ) diff --git a/src/coreComponents/fieldSpecification/SourceFluxStatistics.cpp b/src/coreComponents/fieldSpecification/SourceFluxStatistics.cpp new file mode 100644 index 00000000000..273f2945a24 --- /dev/null +++ b/src/coreComponents/fieldSpecification/SourceFluxStatistics.cpp @@ -0,0 +1,362 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2018-2020 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2020 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2020 TotalEnergies + * Copyright (c) 2019- GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file SourceFluxStatistics.cpp + */ + +#include "SourceFluxStatistics.hpp" + +#include "SourceFluxBoundaryCondition.hpp" +#include "FieldSpecificationManager.hpp" +#include "LvArray/src/tensorOps.hpp" + +namespace geos +{ +using namespace dataRepository; + +SourceFluxStatsAggregator::SourceFluxStatsAggregator( const string & name, + Group * const parent ): + Base( name, parent ) +{ + getWrapper< integer >( Group::viewKeyStruct::logLevelString() ). + appendDescription( GEOS_FMT( "\n- Log Level 1 outputs the sum of all {0}(s) produced rate & mass,\n" + "- Log Level 2 details values for each {0},\n" + "- Log Level 3 details values for each region.", + SourceFluxBoundaryCondition::catalogName() ) ); + + registerWrapper( viewKeyStruct::fluxNamesString().data(), &m_fluxNames ). + setRTTypeName( rtTypes::CustomTypes::groupNameRefArray ). + setInputFlag( InputFlags::OPTIONAL ). + setSizedFromParent( 0 ). + setDefaultValue( "*" ). + setDescription( GEOS_FMT( "Name(s) array of the {0}(s) for which we want the statistics. " + "Use \"*\" to target all {0}.", + SourceFluxBoundaryCondition::catalogName() ) ); +} + +void SourceFluxStatsAggregator::postProcessInput() +{ + Base::postProcessInput(); + + FieldSpecificationManager & fsManager = FieldSpecificationManager::getInstance(); + if( m_fluxNames.size() == 1 && m_fluxNames[0] == "*" ) + { + m_fluxNames.clear(); + fsManager.forSubGroups< SourceFluxBoundaryCondition >( [&]( SourceFluxBoundaryCondition & sourceFlux ) + { + m_fluxNames.emplace_back( string( sourceFlux.getName() ) ); + } ); + GEOS_WARNING_IF( m_fluxNames.empty(), + GEOS_FMT( "{}: No {} was found in {}.", + getDataContext(), SourceFluxBoundaryCondition::catalogName(), + fsManager.getDataContext() ) ); + } + else + { + for( string const & fluxName : m_fluxNames ) + { + GEOS_ERROR_IF( !fsManager.hasGroup< SourceFluxBoundaryCondition >( fluxName ), + GEOS_FMT( "{}: No {} named {} was found in {}.", + getDataContext(), SourceFluxBoundaryCondition::catalogName(), + fluxName, fsManager.getDataContext() ) ); + } + } +} + +Wrapper< SourceFluxStatsAggregator::WrappedStats > & +SourceFluxStatsAggregator::registerWrappedStats( Group & group, + string_view fluxName, + string_view elementSetName ) +{ + string const wrapperName = getStatWrapperName( fluxName ); + Wrapper< WrappedStats > & statsWrapper = group.registerWrapper< WrappedStats >( wrapperName ); + statsWrapper.setRestartFlags( RestartFlags::NO_WRITE ); + statsWrapper.reference().setTarget( getName(), fluxName ); + + writeStatsToCSV( elementSetName, statsWrapper.reference(), true ); + + return statsWrapper; +} +void SourceFluxStatsAggregator::registerDataOnMesh( Group & meshBodies ) +{ + if( m_solver == nullptr ) + { + return; + } + + m_solver->forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, + MeshLevel & mesh, + arrayView1d< string const > const & ) + { + registerWrappedStats( mesh, viewKeyStruct::fluxSetWrapperString(), viewKeyStruct::allRegionWrapperString() ); + + for( string const & fluxName : m_fluxNames ) + { + registerWrappedStats( mesh, fluxName, viewKeyStruct::allRegionWrapperString() ); + + mesh.getElemManager().forElementRegions( [&]( ElementRegionBase & region ) + { + Wrapper< WrappedStats > & regionStatsWrapper = + registerWrappedStats( region, fluxName, region.getName() ); + region.excludeWrappersFromPacking( { regionStatsWrapper.getName() } ); + + region.forElementSubRegions( [&]( ElementSubRegionBase & subRegion ) + { + Wrapper< WrappedStats > & subRegionStatsWrapper = + registerWrappedStats( subRegion, fluxName, subRegion.getName() ); + subRegion.excludeWrappersFromPacking( { subRegionStatsWrapper.getName() } ); + } ); + } ); + } + } ); +} + +void SourceFluxStatsAggregator::writeStatsToLog( integer minLogLevel, + string_view elementSetName, + WrappedStats const & wrappedStats ) +{ + if( getLogLevel() >= minLogLevel && logger::internal::rank == 0 ) + { + GEOS_LOG_RANK( GEOS_FMT( "{} {} (of {}, in {}): Producing on {} elements", + catalogName(), getName(), wrappedStats.getFluxName(), elementSetName, + wrappedStats.stats().m_elementCount ) ); + + // we want to format differently if we have got multiple phases or not + string_view massUnit = units::getSymbol( m_solver->getMassUnit() ); + if( wrappedStats.stats().m_producedMass.size() == 1 ) + { + GEOS_LOG_RANK( GEOS_FMT( "{} {} (of {}, in {}): Produced mass = {} {}", + catalogName(), getName(), wrappedStats.getFluxName(), elementSetName, + wrappedStats.stats().m_producedMass[0], massUnit ) ); + GEOS_LOG_RANK( GEOS_FMT( "{} {} (of {}, in {}): Production rate = {} {}/s", + catalogName(), getName(), wrappedStats.getFluxName(), elementSetName, + wrappedStats.stats().m_productionRate[0], massUnit ) ); + } + else + { + GEOS_LOG_RANK( GEOS_FMT( "{} {} (of {}, in {}): Produced mass = {} {}", + catalogName(), getName(), wrappedStats.getFluxName(), elementSetName, + wrappedStats.stats().m_producedMass, massUnit ) ); + GEOS_LOG_RANK( GEOS_FMT( "{} {} (of {}, in {}): Production rate = {} {}/s", + catalogName(), getName(), wrappedStats.getFluxName(), elementSetName, + wrappedStats.stats().m_productionRate, massUnit ) ); + } + } +} + +void SourceFluxStatsAggregator::writeStatsToCSV( string_view elementSetName, WrappedStats const & stats, + bool writeHeader ) +{ + if( m_writeCSV > 0 && MpiWrapper::commRank() == 0 ) + { + string const fileName = GEOS_FMT( "{}/{}_{}_{}.csv", + m_outputDir, + stats.getAggregatorName(), stats.getFluxName(), elementSetName ); + std::ofstream outputFile( fileName, + writeHeader ? std::ios_base::out : std::ios_base::app ); + if( writeHeader ) + { + outputFile << GEOS_FMT( "Time [s],Element Count,Producted Mass [{0}],Production Rate [{0}/s]", + units::getSymbol( m_solver->getMassUnit() ) ) << std::endl; + } + else + { + outputFile << GEOS_FMT( "{},{},{},{}", + stats.getStatsPeriodStart(), stats.stats().m_elementCount, + stats.stats().m_producedMass, stats.stats().m_productionRate ) << std::endl; + } + outputFile.close(); + } +} + +bool SourceFluxStatsAggregator::execute( real64 const GEOS_UNUSED_PARAM( time_n ), + real64 const GEOS_UNUSED_PARAM( dt ), + integer const GEOS_UNUSED_PARAM( cycleNumber ), + integer const GEOS_UNUSED_PARAM( eventCounter ), + real64 const GEOS_UNUSED_PARAM( eventProgress ), + DomainPartition & domain ) +{ + forMeshLevelStatsWrapper( domain, + [&] ( MeshLevel & meshLevel, WrappedStats & meshLevelStats ) + { + meshLevelStats.stats() = StatData(); + + forAllFluxStatsWrappers( meshLevel, + [&] ( MeshLevel &, WrappedStats & fluxStats ) + { + fluxStats.stats() = StatData(); + + forAllRegionStatsWrappers( meshLevel, fluxStats.getFluxName(), + [&] ( ElementRegionBase & region, WrappedStats & regionStats ) + { + regionStats.stats() = StatData(); + + forAllSubRegionStatsWrappers( region, regionStats.getFluxName(), + [&] ( ElementSubRegionBase &, WrappedStats & subRegionStats ) + { + subRegionStats.finalizePeriod(); + regionStats.stats().combine( subRegionStats.stats() ); + } ); + + fluxStats.stats().combine( regionStats.stats() ); + writeStatsToLog( 3, region.getName(), regionStats ); + writeStatsToCSV( region.getName(), regionStats, false ); + } ); + + meshLevelStats.stats().combine( fluxStats.stats() ); + writeStatsToLog( 2, viewKeyStruct::allRegionWrapperString(), fluxStats ); + writeStatsToCSV( viewKeyStruct::allRegionWrapperString(), fluxStats, false ); + } ); + + writeStatsToLog( 1, viewKeyStruct::allRegionWrapperString(), meshLevelStats ); + writeStatsToCSV( viewKeyStruct::allRegionWrapperString(), meshLevelStats, false ); + } ); + + return false; +} + + + +void SourceFluxStatsAggregator::StatData::allocate( integer phaseCount ) +{ + if( m_producedMass.size() != phaseCount ) + { + m_producedMass.resize( phaseCount ); + m_productionRate.resize( phaseCount ); + } +} +void SourceFluxStatsAggregator::StatData::reset() +{ + for( int ip = 0; ip < getPhaseCount(); ++ip ) + { + m_producedMass[ip] = 0.0; + m_productionRate[ip] = 0.0; + } + m_elementCount = 0; +} +void SourceFluxStatsAggregator::StatData::combine( StatData const & other ) +{ + allocate( other.getPhaseCount() ); + + for( int ip = 0; ip < other.getPhaseCount(); ++ip ) + { + m_producedMass[ip] += other.m_producedMass[ip]; + m_productionRate[ip] += other.m_productionRate[ip]; + } + m_elementCount += other.m_elementCount; +} +void SourceFluxStatsAggregator::StatData::mpiReduce() +{ + for( int ip = 0; ip < getPhaseCount(); ++ip ) + { + m_producedMass[ip] = MpiWrapper::sum( m_producedMass[ip] ); + m_productionRate[ip] = MpiWrapper::sum( m_productionRate[ip] ); + } + m_elementCount = MpiWrapper::sum( m_elementCount ); +} + +void SourceFluxStatsAggregator::WrappedStats::setTarget( string_view aggregatorName, + string_view fluxName ) +{ + m_aggregatorName = aggregatorName; + m_fluxName = fluxName; +} +void SourceFluxStatsAggregator::WrappedStats::gatherTimeStepStats( real64 const currentTime, real64 const dt, + arrayView1d< real64 const > const & producedMass, + integer const elementCount ) +{ + m_periodStats.allocate( producedMass.size() ); + + if( !m_periodStats.m_isGathering ) + { + // if beginning a new period, we must initialize constant values over the period + m_periodStats.m_periodStart = currentTime; + m_periodStats.m_elementCount = elementCount; + m_periodStats.m_isGathering = true; + } + else + { + GEOS_WARNING_IF( currentTime< m_periodStats.m_timeStepStart, GEOS_FMT( "{}: Time seems to have rollback, stats will be wrong.", m_aggregatorName ) ); + if( currentTime > m_periodStats.m_timeStepStart ) + { + // if beginning a new timestep, we must accumulate the stats from previous timesteps (mass & dt) before collecting the new ones + for( int ip = 0; ip < m_periodStats.getPhaseCount(); ++ip ) + { + m_periodStats.m_periodPendingMass[ip] += m_periodStats.m_timeStepMass[ip]; + } + m_periodStats.m_periodPendingDeltaTime += m_periodStats.m_timeStepDeltaTime; + } + } + // current timestep stats to take into account (overriding if not begining a new timestep) + m_periodStats.m_timeStepStart = currentTime; + m_periodStats.m_timeStepDeltaTime = dt; + for( int ip = 0; ip < m_periodStats.getPhaseCount(); ++ip ) + { + m_periodStats.m_timeStepMass = producedMass; + } +} +void SourceFluxStatsAggregator::WrappedStats::finalizePeriod() +{ + // init phase data memory allocation if needed + m_stats.allocate( m_periodStats.getPhaseCount() ); + + // produce the period stats of this rank + m_stats.m_elementCount = m_periodStats.m_elementCount; + m_statsPeriodStart = m_periodStats.m_periodStart; + m_statsPeriodDT = m_periodStats.m_timeStepDeltaTime + m_periodStats.m_periodPendingDeltaTime; + + real64 const timeDivisor = m_statsPeriodDT > 0.0 ? 1.0 / m_statsPeriodDT : 0.0; + for( int ip = 0; ip < m_periodStats.getPhaseCount(); ++ip ) + { + real64 periodMass = m_periodStats.m_timeStepMass[ip] + m_periodStats.m_periodPendingMass[ip]; + m_stats.m_producedMass[ip] = periodMass; + m_stats.m_productionRate[ip] = periodMass * timeDivisor; + } + + // combine period results from all MPI ranks + m_stats.mpiReduce(); + + // start a new timestep + m_periodStats.reset(); +} +void SourceFluxStatsAggregator::WrappedStats::PeriodStats::allocate( integer phaseCount ) +{ + if( m_timeStepMass.size() != phaseCount ) + { + m_timeStepMass.resize( phaseCount ); + m_periodPendingMass.resize( phaseCount ); + } +} +void SourceFluxStatsAggregator::WrappedStats::PeriodStats::reset() +{ + for( int ip = 0; ip < getPhaseCount(); ++ip ) + { + m_timeStepMass[ip] = 0.0; + m_periodPendingMass[ip] = 0.0; + } + m_periodPendingDeltaTime = 0.0; + m_elementCount = 0; + m_timeStepStart = 0.0; + m_timeStepDeltaTime = 0.0; + m_isGathering = false; +} + + +REGISTER_CATALOG_ENTRY( TaskBase, + SourceFluxStatsAggregator, + string const &, + dataRepository::Group * const ) + +} /* namespace geos */ diff --git a/src/coreComponents/fieldSpecification/SourceFluxStatistics.hpp b/src/coreComponents/fieldSpecification/SourceFluxStatistics.hpp new file mode 100644 index 00000000000..e4110c92d36 --- /dev/null +++ b/src/coreComponents/fieldSpecification/SourceFluxStatistics.hpp @@ -0,0 +1,420 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2018-2020 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2020 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2020 TotalEnergies + * Copyright (c) 2019- GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file SourceFluxStatistics.hpp + */ + +#ifndef SRC_CORECOMPONENTS_PHYSICSSOLVERS_FLUIDFLOW_SOURCEFLUXSTATISTICS_HPP_ +#define SRC_CORECOMPONENTS_PHYSICSSOLVERS_FLUIDFLOW_SOURCEFLUXSTATISTICS_HPP_ + +#include "physicsSolvers/FieldStatisticsBase.hpp" +#include "physicsSolvers/fluidFlow/FlowSolverBase.hpp" +#include "mesh/DomainPartition.hpp" + + +namespace geos +{ + +/** + * @class SourceFluxStatsAggregator + * + * Task class allowing for the computation of aggregate statistics of SourceFluxBoundaryCondition + */ +class SourceFluxStatsAggregator final : public FieldStatisticsBase< FlowSolverBase > +{ +public: + + /** + * @brief Aggregated flux statistics data. + */ + struct StatData + { + /// Amount of fluid produced by the flux(es). Negative if injecting. One value for each fluid phase. + /// In kg by default, or in mol if useMass = 0 on the solver. + array1d< real64 > m_producedMass; + /// Flux(es) production rate. Negative if injecting. One value for each fluid phase. + /// In kg/s by default, or in mol/s if useMass = 0 on the solver. + array1d< real64 > m_productionRate; + /// Number of elements in which we are producing / injecting + integer m_elementCount = 0; + + /** + * @brief resize the phase data arrays if needed + * @param phaseCount the count of phases + */ + void allocate( integer phaseCount ); + /** + * @brief reset the stats to 0.0 + */ + void reset(); + /** + * @return the phase count for phasic stats + */ + integer getPhaseCount() const + { return m_producedMass.size(); } + /** + * @brief Aggregate the statistics of the instance with those of another one. + * @param other the other stats structure. + */ + void combine( StatData const & other ); + /** + * @brief Aggregate the statistics of the instance with those from all instances from other MPI ranks. + * Must be called only once per timestep. + */ + void mpiReduce(); + }; + /** + * @brief Class that aggregate statistics of a flux over multiple time-steps for a given + * SourceFluxStatsAggregator and a for a given mesh part (i.e. a subregion, a region...). + */ + class WrappedStats + { +public: + /** + * @brief Set the subjects targeted by the stats. + * @param aggregatorName the name of the targeted SourceFluxStatsAggregator + * @param fluxName the name of the targeted SourceFluxBoundaryCondition + */ + void setTarget( string_view aggregatorName, string_view fluxName ); + + /** + * @brief Set the current time step stats. Accumulate the statistics only if the time is strictly + * after the previous time + dt (override the statistics if the current timestep is cut). + * @param currentTime time of the timestep start since simulation starting + * @param dt time delta of the current timestep + * @param producedMass Amount of produced fluid (see StatData::m_producedMass). + * @param elementCount number of cell elements concerned by this instance + */ + void gatherTimeStepStats( real64 currentTime, real64 dt, + arrayView1d< real64 const > const & producedMass, + integer elementCount ); + + /** + * @brief Finalize the period statistics of each timestep gathering and render data over all mpi ranks. + * The result can be read by the data() assessor. + * @note This method must be synchronously called by all MPI ranks. + */ + void finalizePeriod(); + + /** + * @return the reference to the wrapped stats data collected over the last period (one timestep or more), computed by finalizePeriod() + */ + StatData & stats() + { return m_stats; } + + /** + * @return the reference to the wrapped stats data collected over the last period (one timestep or more), computed by finalizePeriod() + */ + StatData const & stats() const + { return m_stats; } + + /** + * @return the start time of the wrapped stats period (in s) + */ + real64 getStatsPeriodStart() const + { return m_statsPeriodStart; } + + /** + * @return the duration of the wrapped stats period (in s) + */ + real64 getStatsPeriodDuration() const + { return m_statsPeriodDT; } + + /** + * @return the name of the SourceFluxStatsAggregator that want to collect data on this instance. + */ + string_view getAggregatorName() const + { return m_aggregatorName; } + + /** + * @return the name of the SourceFluxBoundaryCondition from which we are collecting data on this instance. + */ + string_view getFluxName() const + { return m_fluxName; } +private: + /// stats data collected over the last period (one timestep or more), computed by finalizePeriod() + StatData m_stats; + /// the start time of the wrapped stats period (in s) + real64 m_statsPeriodStart; + /// the duration of the wrapped stats period (in s) + real64 m_statsPeriodDT; + + /** + * @brief This struct is used to accumulate statistics along one or more timesteps. + */ + struct PeriodStats + { + /// Fluid production during current time-step. Same unit as StatData::productedMass. + array1d< real64 > m_timeStepMass; + /// Fluid production during all previous time-step of the current period. Same unit as StatData::productedMass. + array1d< real64 > m_periodPendingMass; + /// time of when the timestep starts (since the simulation start). + real64 m_timeStepStart = 0.0; + /// time that the current timestep is simulating. + real64 m_timeStepDeltaTime = 0.0; + /// start time of the current period. + real64 m_periodStart = 0.0; + /// delta time from all previous time-step of the current period. + real64 m_periodPendingDeltaTime = 0.0; + /// number of cell elements targeted by this instance + integer m_elementCount = 0; + /// Did the period statistics gathering started ? + bool m_isGathering = false; + + /** + * @brief resize the phase data arrays if needed + * @param phaseCount the count of phases + */ + void allocate( integer phaseCount ); + /** + * @brief reset the stats to 0.0 + */ + void reset(); + /** + * @return the phase count for phasic stats + */ + integer getPhaseCount() const + { return m_timeStepMass.size(); } + } m_periodStats; + + /// Name of the SourceFluxStatsAggregator that want to collect data on this instance. + string m_aggregatorName; + /// Name of the SourceFluxBoundaryCondition from which we are collecting data on this instance. + string m_fluxName; + }; + + + /** + * @brief Constructor for the statistics class + * @param[in] name the name of the task coming from the xml + * @param[in] parent the parent group of the task + */ + SourceFluxStatsAggregator( const string & name, + Group * const parent ); + + /** + * @return The catalog name + */ + static string catalogName() { return "SourceFluxStatistics"; } + + /** + * @copydoc ExecutableGroup::execute() + */ + virtual bool execute( real64 const time_n, + real64 const dt, + integer const cycleNumber, + integer const eventCounter, + real64 const eventProgress, + DomainPartition & domain ) override; + + /** + * @brief Apply a functor to WrappedStats that combines all stats for each target solver + * discretization mesh levels. + * @param domain the domain for which we want the statistics + * @param lambda the functor that will be called for each WrappedStats. Takes in parameter the MeshLevel + * reference and the reference to the WrappedStats that combines all stats for the instance. + * @tparam LAMBDA the type of lambda function to call in the function + * @note To be retrieved, the WrappedStats structs must be registered on the container during the + * registerDataOnMesh() call. + */ + template< typename LAMBDA > + void forMeshLevelStatsWrapper( DomainPartition & domain, LAMBDA && lambda ); + /** + * @brief Apply a functor to each WrappedStats that combines the stats over all region for a flux. + * @param meshLevel the mesh level for which we want the statistics. + * @param lambda the functor that will be called for each WrappedStats. Takes in parameter the MeshLevel + * reference and the reference to one of the flux WrappedStats. + * @tparam LAMBDA the type of lambda function to call in the function + * @note To be retrieved, the WrappedStats structs must be registered on the container during the + * registerDataOnMesh() call. + */ + template< typename LAMBDA > + void forAllFluxStatsWrappers( MeshLevel & meshLevel, LAMBDA && lambda ); + /** + * @brief Apply a functor to all simulated region WrappedStats (of the given MeshLevel) that target a + * given flux. + * @param meshLevel the mesh level we want to loop over all its regions. + * @param fluxName the name of the flux from which we want the statistics. + * @param lambda the functor that will be called for each WrappedStats. Takes in parameter the + * ElementRegionBase reference and the reference of its WrappedStats. + * @tparam LAMBDA the type of lambda function to call in the function + * @note To be retrieved, the WrappedStats structs must be registered on the container during the + * registerDataOnMesh() call. + */ + template< typename LAMBDA > + void forAllRegionStatsWrappers( MeshLevel & meshLevel, string_view fluxName, LAMBDA && lambda ); + /** + * @brief Apply a functor to all subregion WrappedStats (of the given region) that target a given flux. + * @param region the region from which we want to execute the lambda for each of its sub-region. + * @param fluxName the name of the flux from which we want the statistics. + * @param lambda the functor that will be called for each WrappedStats. Takes in parameter the + * ElementSubRegionBase reference and the reference of its WrappedStats. + * @tparam LAMBDA the type of lambda function to call in the function + * @note To be retrieved, the WrappedStats structs must be registered on the container during the + * registerDataOnMesh() call. + */ + template< typename LAMBDA > + void forAllSubRegionStatsWrappers( ElementRegionBase & region, string_view fluxName, LAMBDA && lambda ); + + /** + * @brief Apply a functor to all WrappedStats of the given group that target a given flux (and + * potentially multiple SourceFluxStatsAggregator). + * The functor takes in parameter the reference to the currently processed WrappedStats. + * @note To be retrieved, the WrappedStats structs must be registered on the container during the + * registerDataOnMesh() call. + * @tparam LAMBDA the type of lambda function to call in the function + * @param container the container from which we want the statistics. + * @param fluxName the name of the flux from which we want the statistics. + * @param lambda the functor that will be called for each WrappedStats. Takes in parameter the + * reference to the currently processed WrappedStats. + */ + template< typename LAMBDA > + static void forAllFluxStatWrappers( Group & container, string_view fluxName, LAMBDA && lambda ); + + /** + * @return a string used to name the wrapper that is added to each region that is simulated by the solver. + * The string is unique within the region for the SourceFluxBoundaryCondition and the SourceFluxStatsAggregator. + * @param fluxName The name of the flux. For the mesh-level global stats, fluxSetWrapperString() can be used. + */ + inline string getStatWrapperName( string_view fluxName ) const; + + + /** + * @brief View keys + */ + struct viewKeyStruct + { + /// @return The key for setName + constexpr inline static string_view fluxNamesString() { return "fluxNames"; } + /// @return The key for statistics wrapper name that targets all region + constexpr inline static string_view allRegionWrapperString() { return "all_regions"; } + /// @return The key for statistics wrapper name that targets all fluxes of the set + constexpr inline static string_view fluxSetWrapperString() { return "flux_set"; } + }; + + +private: + using Base = FieldStatisticsBase< FlowSolverBase >; + + /// the names of the SourceFlux(s) for which we want the statistics + string_array m_fluxNames; + + /** + * @copydoc Group::registerDataOnMesh(Group &) + */ + void registerDataOnMesh( Group & meshBodies ) override; + + /** + * @copydoc Group::postProcessInput() + */ + void postProcessInput() override; + + dataRepository::Wrapper< WrappedStats > & registerWrappedStats( Group & group, + string_view fluxName, + string_view elementSetName ); + + /** + * @brief If requested, output in the log and CSV the given statistics. + * @param minLogLevel the min log level to output any line. + * @param elementSetName the region / sub-subregion name concerned by the statistics. + * @param stats the statistics that must be output in the log. + */ + void writeStatsToLog( integer minLogLevel, string_view elementSetName, WrappedStats const & stats ); + /** + * @brief If CSV is enabled, create or append a new CSV file. + * @param elementSetName the region / sub-subregion name concerned by the statistics. + * @param stats the statistics that must be output in the log. + * @param writeHeader If true, create the CSV with the header. If false, append it with the statistics. + */ + void writeStatsToCSV( string_view elementSetName, WrappedStats const & stats, bool writeHeader ); + +}; + + +template< typename LAMBDA > +void SourceFluxStatsAggregator::forAllFluxStatWrappers( Group & container, + string_view fluxName, + LAMBDA && lambda ) +{ + container.forWrappers< WrappedStats >( [&]( dataRepository::Wrapper< WrappedStats > & statsWrapper ) + { + if( statsWrapper.referenceAsView().getFluxName() == fluxName ) + { + lambda( statsWrapper.reference() ); + } + } ); +} + +template< typename LAMBDA > +void SourceFluxStatsAggregator::forMeshLevelStatsWrapper( DomainPartition & domain, + LAMBDA && lambda ) +{ + m_solver->forDiscretizationOnMeshTargets( domain.getMeshBodies(), + [&] ( string const &, + MeshLevel & meshLevel, + arrayView1d< string const > const & ) + { + string const wrapperName = getStatWrapperName( viewKeyStruct::fluxSetWrapperString() ); + WrappedStats & stats = meshLevel.getWrapper< WrappedStats >( wrapperName ).reference(); + + lambda( meshLevel, stats ); + } ); +} +template< typename LAMBDA > +void SourceFluxStatsAggregator::forAllFluxStatsWrappers( MeshLevel & meshLevel, + LAMBDA && lambda ) +{ + for( string const & fluxName : m_fluxNames ) + { + string const wrapperName = getStatWrapperName( fluxName ); + WrappedStats & stats = meshLevel.getWrapper< WrappedStats >( wrapperName ).reference(); + + lambda( meshLevel, stats ); + } +} +template< typename LAMBDA > +void SourceFluxStatsAggregator::forAllRegionStatsWrappers( MeshLevel & meshLevel, + string_view fluxName, + LAMBDA && lambda ) +{ + string const wrapperName = getStatWrapperName( fluxName ); + meshLevel.getElemManager().forElementRegions( [&]( ElementRegionBase & region ) + { + WrappedStats & stats = region.getWrapper< WrappedStats >( wrapperName ).reference(); + + lambda( region, stats ); + } ); +} +template< typename LAMBDA > +void SourceFluxStatsAggregator::forAllSubRegionStatsWrappers( ElementRegionBase & region, + string_view fluxName, + LAMBDA && lambda ) +{ + string const wrapperName = getStatWrapperName( fluxName ); + region.forElementSubRegions( [&]( ElementSubRegionBase & subRegion ) + { + WrappedStats & stats = subRegion.getWrapper< WrappedStats >( wrapperName ).reference(); + + lambda( subRegion, stats ); + } ); +} + +inline string SourceFluxStatsAggregator::getStatWrapperName( string_view fluxName ) const +{ return GEOS_FMT( "{}_region_stats_for_{}", fluxName, getName() ); } + + +} /* namespace geos */ + +#endif /* SRC_CORECOMPONENTS_PHYSICSSOLVERS_FLUIDFLOW_SOURCEFLUXSTATISTICS_HPP_ */ diff --git a/src/coreComponents/mesh/ElementRegionBase.hpp b/src/coreComponents/mesh/ElementRegionBase.hpp index 049173a74e6..15a46e407b9 100644 --- a/src/coreComponents/mesh/ElementRegionBase.hpp +++ b/src/coreComponents/mesh/ElementRegionBase.hpp @@ -214,6 +214,20 @@ class ElementRegionBase : public ObjectManagerBase template< typename CONSTITUTIVE_TYPE > string_array getConstitutiveNames() const; + /** + * @return the parent region of the given sub-region. + * @param subRegion the sub-region that we want the parent. + */ + static ElementRegionBase & getParentRegion( ElementSubRegionBase & subRegion ) + { return dynamicCast< ElementRegionBase & >( subRegion.getParent().getParent() ); } + + /** + * @return the parent region of the given sub-region. + * @param subRegion the sub-region that we want the parent. + */ + static ElementRegionBase const & getParentRegion( ElementSubRegionBase const & subRegion ) + { return dynamicCast< ElementRegionBase const & >( subRegion.getParent().getParent() ); } + ///@} diff --git a/src/coreComponents/physicsSolvers/SolverBase.hpp b/src/coreComponents/physicsSolvers/SolverBase.hpp index f9e07f2dc6b..7c47c07639d 100644 --- a/src/coreComponents/physicsSolvers/SolverBase.hpp +++ b/src/coreComponents/physicsSolvers/SolverBase.hpp @@ -734,6 +734,7 @@ class SolverBase : public ExecutableGroup virtual bool registerCallback( void * func, const std::type_info & funcType ) final override; SolverStatistics & getSolverStatistics() { return m_solverStatistics; } + SolverStatistics const & getSolverStatistics() const { return m_solverStatistics; } /** * @brief Return PySolver type. diff --git a/src/coreComponents/physicsSolvers/SolverStatistics.hpp b/src/coreComponents/physicsSolvers/SolverStatistics.hpp index 7ec3b7cc532..964cbc078ca 100644 --- a/src/coreComponents/physicsSolvers/SolverStatistics.hpp +++ b/src/coreComponents/physicsSolvers/SolverStatistics.hpp @@ -79,6 +79,54 @@ class SolverStatistics : public dataRepository::Group */ void outputStatistics() const; + /** + * @return Number of time steps + */ + integer getNumTimeSteps() const + { return m_numTimeSteps; } + + /** + * @return Number of time step cuts + */ + integer getNumTimeStepCuts() const + { return m_numTimeStepCuts; } + + /** + * @return Cumulative number of successful outer loop iterations + */ + integer getNumSuccessfulOuterLoopIterations() const + { return m_numSuccessfulOuterLoopIterations; } + + /** + * @return Cumulative number of successful nonlinear iterations + */ + integer getNumSuccessfulNonlinearIterations() const + { return m_numSuccessfulNonlinearIterations; } + + /** + * @return Cumulative number of successful linear iterations + */ + integer getNumSuccessfulLinearIterations() const + { return m_numSuccessfulLinearIterations; } + + /** + * @return Cumulative number of discarded outer loop iterations + */ + integer getNumDiscardedOuterLoopIterations() const + { return m_numDiscardedOuterLoopIterations; } + + /** + * @return Cumulative number of discarded nonlinear iterations + */ + integer getNumDiscardedNonlinearIterations() const + { return m_numDiscardedNonlinearIterations; } + + /** + * @return Cumulative number of discarded linear iterations + */ + integer getNumDiscardedLinearIterations() const + { return m_numDiscardedLinearIterations; } + private: /** diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.cpp index 6129754a029..652682cca12 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.cpp @@ -35,6 +35,7 @@ #include "fieldSpecification/AquiferBoundaryCondition.hpp" #include "fieldSpecification/EquilibriumInitialCondition.hpp" #include "fieldSpecification/SourceFluxBoundaryCondition.hpp" +#include "fieldSpecification/SourceFluxStatistics.hpp" #include "mesh/DomainPartition.hpp" #include "mesh/mpiCommunications/CommunicationTools.hpp" #include "physicsSolvers/fluidFlow/CompositionalMultiphaseBaseFields.hpp" @@ -73,10 +74,13 @@ CompositionalMultiphaseBase::CompositionalMultiphaseBase( const string & name, setInputFlag( InputFlags::REQUIRED ). setDescription( "Temperature" ); //END_SPHINX_INCLUDE_00 + + // TODO: add more description on how useMass can alter the simulation (convergence issues?). How does it interact with wells useMass? this->registerWrapper( viewKeyStruct::useMassFlagString(), &m_useMass ). setApplyDefaultValue( 0 ). setInputFlag( InputFlags::OPTIONAL ). - setDescription( "Use mass formulation instead of molar" ); + setDescription( GEOS_FMT( "Use mass formulation instead of molar. Warning : Affects {} rates units.", + SourceFluxBoundaryCondition::catalogName() ) ); this->registerWrapper( viewKeyStruct::solutionChangeScalingFactorString(), &m_solutionChangeScalingFactor ). setSizedFromParent( 0 ). @@ -1552,6 +1556,8 @@ void CompositionalMultiphaseBase::applySourceFluxBC( real64 const time, arrayView1d< real64 > rhsContributionArrayView = rhsContributionArray.toView(); localIndex const rankOffset = dofManager.rankOffset(); + RAJA::ReduceSum< parallelDeviceReduce, real64 > massProd( 0.0 ); + // note that the dofArray will not be used after this step (simpler to use dofNumber instead) fs.computeRhsContribution< FieldSpecificationAdd, parallelDevicePolicy<> >( targetSet.toViewConst(), @@ -1585,7 +1591,8 @@ void CompositionalMultiphaseBase::applySourceFluxBC( real64 const time, useTotalMassEquation, dofNumber, rhsContributionArrayView, - localRhs] GEOS_HOST_DEVICE ( localIndex const a ) + localRhs, + massProd] GEOS_HOST_DEVICE ( localIndex const a ) { // we need to filter out ghosts here, because targetSet may contain them localIndex const ei = targetSet[a]; @@ -1594,27 +1601,34 @@ void CompositionalMultiphaseBase::applySourceFluxBC( real64 const time, return; } + real64 const rhsValue = rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the sizeScalingFactor here! + massProd += rhsValue; if( useTotalMassEquation > 0 ) { // for all "fluid components", we add the value to the total mass balance equation globalIndex const totalMassBalanceRow = dofNumber[ei] - rankOffset; - localRhs[totalMassBalanceRow] += rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the - // sizeScalingFactor - // here + localRhs[totalMassBalanceRow] += rhsValue; if( fluidComponentId < numFluidComponents - 1 ) { globalIndex const compMassBalanceRow = totalMassBalanceRow + fluidComponentId + 1; // component mass bal equations are shifted - localRhs[compMassBalanceRow] += rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the - // sizeScalingFactor here + localRhs[compMassBalanceRow] += rhsValue; } } else { globalIndex const compMassBalanceRow = dofNumber[ei] - rankOffset + fluidComponentId; - localRhs[compMassBalanceRow] += rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the - // sizeScalingFactor here + localRhs[compMassBalanceRow] += rhsValue; } } ); + + SourceFluxStatsAggregator::forAllFluxStatWrappers( subRegion, fs.getName(), + [&]( SourceFluxStatsAggregator::WrappedStats & wrapper ) + { + // set the new sub-region statistics for this timestep + array1d< real64 > massProdArr{ m_numComponents }; + massProdArr[fluidComponentId] = massProd.get(); + wrapper.gatherTimeStepStats( time, dt, massProdArr.toViewConst(), targetSet.size() ); + } ); } ); } ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.hpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.hpp index 870213176a7..0708804af11 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.hpp @@ -183,6 +183,12 @@ class CompositionalMultiphaseBase : public FlowSolverBase */ string referenceFluidModelName() const { return m_referenceFluidModelName; } + /** + * @return The unit in which we evaluate the amount of fluid per element (Mass or Mole, depending on useMass). + */ + virtual units::Unit getMassUnit() const override + { return m_useMass ? units::Unit::Mass : units::Unit::Mole; } + /** * @brief assembles the accumulation and volume balance terms for all cells * @param time_n previous time value diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseStatistics.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseStatistics.cpp index 137a760069d..ab5932ce13a 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseStatistics.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseStatistics.cpp @@ -115,8 +115,7 @@ void CompositionalMultiphaseStatistics::registerDataOnMesh( Group & meshBodies ) if( m_writeCSV > 0 && MpiWrapper::commRank() == 0 ) { std::ofstream outputFile( m_outputDir + "/" + regionNames[i] + ".csv" ); - integer const useMass = m_solver->getReference< integer >( CompositionalMultiphaseBase::viewKeyStruct::useMassFlagString() ); - string const massUnit = useMass ? "kg" : "mol"; + string_view massUnit = units::getSymbol( m_solver->getMassUnit() ); outputFile << "Time [s],Min pressure [Pa],Average pressure [Pa],Max pressure [Pa],Min delta pressure [Pa],Max delta pressure [Pa]," << "Min temperature [Pa],Average temperature [Pa],Max temperature [Pa],Total dynamic pore volume [rm^3]"; @@ -305,7 +304,7 @@ void CompositionalMultiphaseStatistics::computeRegionStatistics( real64 const ti subRegionImmobilePhaseMass.toView(), subRegionComponentMass.toView() ); - ElementRegionBase & region = elemManager.getRegion( subRegion.getParent().getParent().getName() ); + ElementRegionBase & region = elemManager.getRegion( ElementRegionBase::getParentRegion( subRegion ).getName() ); RegionStatistics & regionStatistics = region.getReference< RegionStatistics >( viewKeyStruct::regionStatisticsString() ); regionStatistics.averagePressure += subRegionAvgPresNumerator; @@ -405,8 +404,7 @@ void CompositionalMultiphaseStatistics::computeRegionStatistics( real64 const ti mobilePhaseMass[ip] = regionStatistics.phaseMass[ip] - regionStatistics.immobilePhaseMass[ip]; } - integer const useMass = m_solver->getReference< integer >( CompositionalMultiphaseBase::viewKeyStruct::useMassFlagString() ); - string const massUnit = useMass ? "kg" : "mol"; + string_view massUnit = units::getSymbol( m_solver->getMassUnit() ); GEOS_LOG_LEVEL_RANK_0( 1, GEOS_FMT( "{}, {} (time {} s): Pressure (min, average, max): {}, {}, {} Pa", getName(), regionNames[i], time, regionStatistics.minPressure, regionStatistics.averagePressure, regionStatistics.maxPressure ) ); diff --git a/src/coreComponents/physicsSolvers/fluidFlow/FlowSolverBase.hpp b/src/coreComponents/physicsSolvers/fluidFlow/FlowSolverBase.hpp index 3de1b21e7ea..a7d0f163085 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/FlowSolverBase.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/FlowSolverBase.hpp @@ -20,6 +20,7 @@ #define GEOS_PHYSICSSOLVERS_FINITEVOLUME_FLOWSOLVERBASE_HPP_ #include "physicsSolvers/SolverBase.hpp" +#include "common/Units.hpp" namespace geos { @@ -135,6 +136,12 @@ class FlowSolverBase : public SolverBase integer & isThermal() { return m_isThermal; } + /** + * @return The unit in which we evaluate the amount of fluid per element (Mass or Mole). + */ + virtual units::Unit getMassUnit() const + { return units::Unit::Mass; } + /** * @brief Function to activate the flag allowing negative pressure */ diff --git a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp index 85f0354c53a..a08168784b6 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp @@ -23,6 +23,7 @@ #include "discretizationMethods/NumericalMethodsManager.hpp" #include "fieldSpecification/FieldSpecificationManager.hpp" #include "fieldSpecification/SourceFluxBoundaryCondition.hpp" +#include "fieldSpecification/SourceFluxStatistics.hpp" #include "finiteVolume/FiniteVolumeManager.hpp" #include "finiteVolume/FluxApproximationBase.hpp" #include "mesh/DomainPartition.hpp" @@ -861,6 +862,8 @@ void ReactiveCompositionalMultiphaseOBL::applySourceFluxBC( real64 const time, arrayView1d< real64 > rhsContributionArrayView = rhsContributionArray.toView(); localIndex const rankOffset = dofManager.rankOffset(); + RAJA::ReduceSum< parallelDeviceReduce, real64 > massProd( 0.0 ); + // note that the dofArray will not be used after this step (simpler to use dofNumber instead) fs.computeRhsContribution< FieldSpecificationAdd, parallelDevicePolicy<> >( targetSet.toViewConst(), @@ -890,7 +893,8 @@ void ReactiveCompositionalMultiphaseOBL::applySourceFluxBC( real64 const time, fluidComponentId, dofNumber, rhsContributionArrayView, - localRhs] GEOS_HOST_DEVICE ( localIndex const a ) + localRhs, + massProd] GEOS_HOST_DEVICE ( localIndex const a ) { // we need to filter out ghosts here, because targetSet may contain them localIndex const ei = targetSet[a]; @@ -901,7 +905,18 @@ void ReactiveCompositionalMultiphaseOBL::applySourceFluxBC( real64 const time, // for all "fluid components", we add the value to the component mass balance equation globalIndex const compMassBalanceRow = dofNumber[ei] - rankOffset + fluidComponentId; - localRhs[compMassBalanceRow] += rhsContributionArrayView[a] / sizeScalingFactor; + real64 const rhsValue = rhsContributionArrayView[a] / sizeScalingFactor; + localRhs[compMassBalanceRow] += rhsValue; + massProd += rhsValue; + } ); + + SourceFluxStatsAggregator::forAllFluxStatWrappers( subRegion, fs.getName(), + [&]( SourceFluxStatsAggregator::WrappedStats & wrapper ) + { + // set the new sub-region statistics for this timestep + array1d< real64 > massProdArr{ m_numComponents }; + massProdArr[fluidComponentId] = massProd.get(); + wrapper.gatherTimeStepStats( time, dt, massProdArr.toViewConst(), targetSet.size() ); } ); } ); } ); diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp index 46066692beb..c342cc8a592 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp @@ -31,6 +31,7 @@ #include "fieldSpecification/EquilibriumInitialCondition.hpp" #include "fieldSpecification/FieldSpecificationManager.hpp" #include "fieldSpecification/SourceFluxBoundaryCondition.hpp" +#include "fieldSpecification/SourceFluxStatistics.hpp" #include "finiteVolume/FiniteVolumeManager.hpp" #include "functions/TableFunction.hpp" #include "mainInterface/ProblemManager.hpp" @@ -57,6 +58,13 @@ SinglePhaseBase::SinglePhaseBase( const string & name, setApplyDefaultValue( 0.0 ). setInputFlag( InputFlags::OPTIONAL ). setDescription( "Temperature" ); + + this->getWrapper< integer >( string( viewKeyStruct::isThermalString() ) ). + appendDescription( GEOS_FMT( "\nSourceFluxes application if {} is enabled :\n" + "- negative value (injection): the mass balance equation is modified to considered the additional source term,\n" + "- positive value (production): both the mass balance and the energy balance equations are modified to considered the additional source term.\n" + "For the energy balance equation, the mass flux is multipied by the enthalpy in the cell from which the fluid is being produced.", + viewKeyStruct::isThermalString() ) ); } @@ -877,20 +885,12 @@ void SinglePhaseBase::applyBoundaryConditions( real64 time_n, namespace { -char const bcLogMessage[] = - "SinglePhaseBase {}: at time {}s, " - "the <{}> boundary condition '{}' is applied to the element set '{}' in subRegion '{}'. " - "\nThe scale of this boundary condition is {} and multiplies the value of the provided function (if any). " - "\nThe total number of target elements (including ghost elements) is {}. " - "\nNote that if this number is equal to zero for all subRegions, the boundary condition will not be applied on this element set."; - void applyAndSpecifyFieldValue( real64 const & time_n, real64 const & dt, MeshLevel & mesh, globalIndex const rankOffset, string const dofKey, - bool const isFirstNonlinearIteration, - string const solverName, + bool const, integer const idof, string const fieldKey, string const boundaryFieldKey, @@ -903,19 +903,11 @@ void applyAndSpecifyFieldValue( real64 const & time_n, mesh, fieldKey, [&]( FieldSpecificationBase const & fs, - string const & setName, + string const &, SortedArrayView< localIndex const > const & lset, ElementSubRegionBase & subRegion, string const & ) { - if( fs.getLogLevel() >= 1 && isFirstNonlinearIteration ) - { - globalIndex const numTargetElems = MpiWrapper::sum< globalIndex >( lset.size() ); - GEOS_LOG_RANK_0( GEOS_FMT( bcLogMessage, - solverName, time_n+dt, fs.getCatalogName(), fs.getName(), - setName, subRegion.getName(), fs.getScale(), numTargetElems ) ); - } - // Specify the bc value of the field fs.applyFieldValue< FieldSpecificationEqual, parallelDevicePolicy<> >( lset, @@ -974,12 +966,12 @@ void SinglePhaseBase::applyDirichletBC( real64 const time_n, MeshLevel & mesh, arrayView1d< string const > const & ) { - applyAndSpecifyFieldValue( time_n, dt, mesh, rankOffset, dofKey, isFirstNonlinearIteration, getName(), + applyAndSpecifyFieldValue( time_n, dt, mesh, rankOffset, dofKey, isFirstNonlinearIteration, 0, fields::flow::pressure::key(), fields::flow::bcPressure::key(), localMatrix, localRhs ); if( m_isThermal ) { - applyAndSpecifyFieldValue( time_n, dt, mesh, rankOffset, dofKey, isFirstNonlinearIteration, getName(), + applyAndSpecifyFieldValue( time_n, dt, mesh, rankOffset, dofKey, isFirstNonlinearIteration, 1, fields::flow::temperature::key(), fields::flow::bcTemperature::key(), localMatrix, localRhs ); } @@ -1044,25 +1036,6 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, ElementSubRegionBase & subRegion, string const & ) { - if( fs.getLogLevel() >= 1 && m_nonlinearSolverParameters.m_numNewtonIterations == 0 ) - { - globalIndex const numTargetElems = MpiWrapper::sum< globalIndex >( targetSet.size() ); - GEOS_LOG_RANK_0( GEOS_FMT( bcLogMessage, - getName(), time_n+dt, fs.getCatalogName(), fs.getName(), - setName, subRegion.getName(), fs.getScale(), numTargetElems ) ); - - if( isThermal ) - { - char const msg[] = "SinglePhaseBase {} with isThermal = 1. At time {}s, " - "the <{}> source flux boundary condition '{}' will be applied with the following behavior" - "\n - negative value (injection): the mass balance equation is modified to considered the additional source term" - "\n - positive value (production): both the mass balance and the energy balance equations are modified to considered the additional source term. " \ - "\n For the energy balance equation, the mass flux is multipied by the enthalpy in the cell from which the fluid is being produced."; - GEOS_LOG_RANK_0( GEOS_FMT( msg, - getName(), time_n+dt, fs.getCatalogName(), fs.getName() ) ); - } - } - if( targetSet.size() == 0 ) { return; @@ -1087,6 +1060,8 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, arrayView1d< real64 > rhsContributionArrayView = rhsContributionArray.toView(); localIndex const rankOffset = dofManager.rankOffset(); + RAJA::ReduceSum< parallelDeviceReduce, real64 > massProd( 0.0 ); + // note that the dofArray will not be used after this step (simpler to use dofNumber instead) fs.computeRhsContribution< FieldSpecificationAdd, parallelDevicePolicy<> >( targetSet.toViewConst(), @@ -1126,7 +1101,8 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, dEnthalpy_dPressure, rhsContributionArrayView, localRhs, - localMatrix] GEOS_HOST_DEVICE ( localIndex const a ) + localMatrix, + massProd] GEOS_HOST_DEVICE ( localIndex const a ) { // we need to filter out ghosts here, because targetSet may contain them localIndex const ei = targetSet[a]; @@ -1140,6 +1116,7 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, globalIndex const energyRowIndex = massRowIndex + 1; real64 const rhsValue = rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the sizeScalingFactor here! localRhs[massRowIndex] += rhsValue; + massProd += rhsValue; //add the value to the energey balance equation if the flux is positive (i.e., it's a producer) if( rhsContributionArrayView[a] > 0.0 ) { @@ -1166,7 +1143,8 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, ghostRank, dofNumber, rhsContributionArrayView, - localRhs] GEOS_HOST_DEVICE ( localIndex const a ) + localRhs, + massProd] GEOS_HOST_DEVICE ( localIndex const a ) { // we need to filter out ghosts here, because targetSet may contain them localIndex const ei = targetSet[a]; @@ -1177,9 +1155,20 @@ void SinglePhaseBase::applySourceFluxBC( real64 const time_n, // add the value to the mass balance equation globalIndex const rowIndex = dofNumber[ei] - rankOffset; - localRhs[rowIndex] += rhsContributionArrayView[a] / sizeScalingFactor; // scale the contribution by the sizeScalingFactor here! + real64 const rhsValue = rhsContributionArrayView[a] / sizeScalingFactor; + localRhs[rowIndex] += rhsValue; + massProd += rhsValue; } ); } + + SourceFluxStatsAggregator::forAllFluxStatWrappers( subRegion, fs.getName(), + [&]( SourceFluxStatsAggregator::WrappedStats & wrapper ) + { + // set the new sub-region statistics for this timestep + array1d< real64 > massProdArr{ 1 }; + massProdArr[0] = massProd.get(); + wrapper.gatherTimeStepStats( time_n, dt, massProdArr.toViewConst(), targetSet.size() ); + } ); } ); } ); } diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.cpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.cpp index c74a9e906b1..73b34b0703b 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.cpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.cpp @@ -175,7 +175,7 @@ void SinglePhaseStatistics::computeRegionStatistics( real64 const time, subRegionTotalPoreVol, subRegionTotalMass ); - ElementRegionBase & region = elemManager.getRegion( subRegion.getParent().getParent().getName() ); + ElementRegionBase & region = elemManager.getRegion( ElementRegionBase::getParentRegion( subRegion ).getName() ); RegionStatistics & regionStatistics = region.getReference< RegionStatistics >( viewKeyStruct::regionStatisticsString() ); regionStatistics.averagePressure += subRegionAvgPresNumerator; diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.hpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.hpp index bfd1360fd92..57ac9fc1f64 100644 --- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.hpp +++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseStatistics.hpp @@ -62,9 +62,6 @@ class SinglePhaseStatistics : public FieldStatisticsBase< SinglePhaseBase > /**@}*/ -private: - - using Base = FieldStatisticsBase< SinglePhaseBase >; /** * @struct viewKeyStruct holds char strings and viewKeys for fast lookup @@ -105,6 +102,10 @@ class SinglePhaseStatistics : public FieldStatisticsBase< SinglePhaseBase > real64 totalUncompactedPoreVolume; }; +private: + + using Base = FieldStatisticsBase< SinglePhaseBase >; + /** * @brief Compute some statistics on the reservoir (average field pressure, etc) * @param[in] mesh the mesh level object diff --git a/src/coreComponents/schema/docs/CompositionalMultiphaseFVM.rst b/src/coreComponents/schema/docs/CompositionalMultiphaseFVM.rst index 1c47dd14fd5..33e2a2e5cb5 100644 --- a/src/coreComponents/schema/docs/CompositionalMultiphaseFVM.rst +++ b/src/coreComponents/schema/docs/CompositionalMultiphaseFVM.rst @@ -36,7 +36,7 @@ targetRelativePressureChangeInTimeStep real64 targetRelativeTemperatureChangeInTimeStep real64 0.2 Target (relative) change in temperature in a time step (expected value between 0 and 1) temperature real64 required Temperature useDBC integer 0 Enable Dissipation-based continuation flux -useMass integer 0 Use mass formulation instead of molar +useMass integer 0 Use mass formulation instead of molar. Warning : Affects SourceFlux rates units. useSimpleAccumulation integer 1 Flag indicating whether simple accumulation form is used useTotalMassEquation integer 1 Flag indicating whether total mass equation is used LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` diff --git a/src/coreComponents/schema/docs/CompositionalMultiphaseHybridFVM.rst b/src/coreComponents/schema/docs/CompositionalMultiphaseHybridFVM.rst index 067a6a06a6a..1627f7eae15 100644 --- a/src/coreComponents/schema/docs/CompositionalMultiphaseHybridFVM.rst +++ b/src/coreComponents/schema/docs/CompositionalMultiphaseHybridFVM.rst @@ -27,7 +27,7 @@ targetRegions groupNameRef_array required Allowable targetRelativePressureChangeInTimeStep real64 0.2 Target (relative) change in pressure in a time step (expected value between 0 and 1) targetRelativeTemperatureChangeInTimeStep real64 0.2 Target (relative) change in temperature in a time step (expected value between 0 and 1) temperature real64 required Temperature -useMass integer 0 Use mass formulation instead of molar +useMass integer 0 Use mass formulation instead of molar. Warning : Affects SourceFlux rates units. useSimpleAccumulation integer 1 Flag indicating whether simple accumulation form is used useTotalMassEquation integer 1 Flag indicating whether total mass equation is used LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` diff --git a/src/coreComponents/schema/docs/SinglePhaseFVM.rst b/src/coreComponents/schema/docs/SinglePhaseFVM.rst index 6000056c728..f6ffdc9cbfb 100644 --- a/src/coreComponents/schema/docs/SinglePhaseFVM.rst +++ b/src/coreComponents/schema/docs/SinglePhaseFVM.rst @@ -1,22 +1,26 @@ -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -Name Type Default Description -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -allowNegativePressure integer 1 Flag indicating if negative pressure is allowed -cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] -discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. -initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. -isThermal integer 0 Flag indicating whether the problem is thermal or not. -logLevel integer 0 Log level -maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration -maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check -maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check -name groupName required A name is required for any non-unique nodes -targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. -temperature real64 0 Temperature -LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` -NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +Name Type Default Description +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +allowNegativePressure integer 1 Flag indicating if negative pressure is allowed +cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] +discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. +initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. +isThermal integer 0 | Flag indicating whether the problem is thermal or not. + | SourceFluxes application if isThermal is enabled : + | - negative value (injection): the mass balance equation is modified to considered the additional source term, + | - positive value (production): both the mass balance and the energy balance equations are modified to considered the additional source term. + | For the energy balance equation, the mass flux is multipied by the enthalpy in the cell from which the fluid is being produced. +logLevel integer 0 Log level +maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration +maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check +maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check +name groupName required A name is required for any non-unique nodes +targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. +temperature real64 0 Temperature +LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` +NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== diff --git a/src/coreComponents/schema/docs/SinglePhaseHybridFVM.rst b/src/coreComponents/schema/docs/SinglePhaseHybridFVM.rst index 6000056c728..f6ffdc9cbfb 100644 --- a/src/coreComponents/schema/docs/SinglePhaseHybridFVM.rst +++ b/src/coreComponents/schema/docs/SinglePhaseHybridFVM.rst @@ -1,22 +1,26 @@ -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -Name Type Default Description -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -allowNegativePressure integer 1 Flag indicating if negative pressure is allowed -cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] -discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. -initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. -isThermal integer 0 Flag indicating whether the problem is thermal or not. -logLevel integer 0 Log level -maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration -maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check -maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check -name groupName required A name is required for any non-unique nodes -targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. -temperature real64 0 Temperature -LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` -NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +Name Type Default Description +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +allowNegativePressure integer 1 Flag indicating if negative pressure is allowed +cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] +discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. +initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. +isThermal integer 0 | Flag indicating whether the problem is thermal or not. + | SourceFluxes application if isThermal is enabled : + | - negative value (injection): the mass balance equation is modified to considered the additional source term, + | - positive value (production): both the mass balance and the energy balance equations are modified to considered the additional source term. + | For the energy balance equation, the mass flux is multipied by the enthalpy in the cell from which the fluid is being produced. +logLevel integer 0 Log level +maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration +maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check +maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check +name groupName required A name is required for any non-unique nodes +targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. +temperature real64 0 Temperature +LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` +NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== diff --git a/src/coreComponents/schema/docs/SinglePhaseProppantFVM.rst b/src/coreComponents/schema/docs/SinglePhaseProppantFVM.rst index 6000056c728..f6ffdc9cbfb 100644 --- a/src/coreComponents/schema/docs/SinglePhaseProppantFVM.rst +++ b/src/coreComponents/schema/docs/SinglePhaseProppantFVM.rst @@ -1,22 +1,26 @@ -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -Name Type Default Description -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== -allowNegativePressure integer 1 Flag indicating if negative pressure is allowed -cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] -discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. -initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. -isThermal integer 0 Flag indicating whether the problem is thermal or not. -logLevel integer 0 Log level -maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration -maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check -maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check -name groupName required A name is required for any non-unique nodes -targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. -temperature real64 0 Temperature -LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` -NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` -============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================== +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +Name Type Default Description +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== +allowNegativePressure integer 1 Flag indicating if negative pressure is allowed +cflFactor real64 0.5 Factor to apply to the `CFL condition `_ when calculating the maximum allowable time step. Values should be in the interval (0,1] +discretization groupNameRef required Name of discretization object (defined in the :ref:`NumericalMethodsManager`) to use for this solver. For instance, if this is a Finite Element Solver, the name of a :ref:`FiniteElement` should be specified. If this is a Finite Volume Method, the name of a :ref:`FiniteVolume` discretization should be specified. +initialDt real64 1e+99 Initial time-step value required by the solver to the event manager. +isThermal integer 0 | Flag indicating whether the problem is thermal or not. + | SourceFluxes application if isThermal is enabled : + | - negative value (injection): the mass balance equation is modified to considered the additional source term, + | - positive value (production): both the mass balance and the energy balance equations are modified to considered the additional source term. + | For the energy balance equation, the mass flux is multipied by the enthalpy in the cell from which the fluid is being produced. +logLevel integer 0 Log level +maxAbsolutePressureChange real64 -1 Maximum (absolute) pressure change in a Newton iteration +maxSequentialPressureChange real64 100000 Maximum (absolute) pressure change in a sequential iteration, used for outer loop convergence check +maxSequentialTemperatureChange real64 0.1 Maximum (absolute) temperature change in a sequential iteration, used for outer loop convergence check +name groupName required A name is required for any non-unique nodes +targetRegions groupNameRef_array required Allowable regions that the solver may be applied to. Note that this does not indicate that the solver will be applied to these regions, only that allocation will occur such that the solver may be applied to these regions. The decision about what regions this solver will beapplied to rests in the EventManager. +temperature real64 0 Temperature +LinearSolverParameters node unique :ref:`XML_LinearSolverParameters` +NonlinearSolverParameters node unique :ref:`XML_NonlinearSolverParameters` +============================== ================== ======== ======================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================== diff --git a/src/coreComponents/schema/docs/SourceFlux.rst b/src/coreComponents/schema/docs/SourceFlux.rst index c227dcbd7d8..660880a3b49 100644 --- a/src/coreComponents/schema/docs/SourceFlux.rst +++ b/src/coreComponents/schema/docs/SourceFlux.rst @@ -1,20 +1,20 @@ -====================== ================== ======== ============================================================================== -Name Type Default Description -====================== ================== ======== ============================================================================== -bcApplicationTableName groupNameRef Name of table that specifies the on/off application of the boundary condition. -beginTime real64 -1e+99 Time at which the boundary condition will start being applied. -component integer -1 Component of field (if tensor) to apply boundary condition to. -direction R1Tensor {0,0,0} Direction to apply boundary condition to. -endTime real64 1e+99 Time at which the boundary condition will stop being applied. -functionName groupNameRef Name of function that specifies variation of the boundary condition. -initialCondition integer 0 Boundary condition is applied as an initial condition. -logLevel integer 0 Log level -name groupName required A name is required for any non-unique nodes -objectPath groupNameRef Path to the target field -scale real64 0 Scale factor for value of the boundary condition. -setNames groupNameRef_array required Name of sets that boundary condition is applied to. -====================== ================== ======== ============================================================================== +====================== ================== ======== ======================================================================================================================================================================================================================================================================================== +Name Type Default Description +====================== ================== ======== ======================================================================================================================================================================================================================================================================================== +bcApplicationTableName groupNameRef Name of table that specifies the on/off application of the boundary condition. +beginTime real64 -1e+99 Time at which the boundary condition will start being applied. +component integer -1 Component of field (if tensor) to apply boundary condition to. +direction R1Tensor {0,0,0} Direction to apply boundary condition to. +endTime real64 1e+99 Time at which the boundary condition will stop being applied. +functionName groupNameRef Name of a function that specifies the variation of the production rate variations of this SourceFlux.Multiplied by scale. If no function is provided, a constant value of 1 is used.The producted fluid rate unit is in kg by default, or in mole if the flow solver has a useMass of 0. +initialCondition integer 0 Boundary condition is applied as an initial condition. +logLevel integer 0 Log level +name groupName required A name is required for any non-unique nodes +objectPath groupNameRef Path to the target field +scale real64 0 Multiplier of the functionName value. If no functionName is provided, this value is used directly. +setNames groupNameRef_array required Name of sets that boundary condition is applied to. +====================== ================== ======== ======================================================================================================================================================================================================================================================================================== diff --git a/src/coreComponents/schema/docs/SourceFluxStatistics.rst b/src/coreComponents/schema/docs/SourceFluxStatistics.rst new file mode 100644 index 00000000000..a5ef4f42a89 --- /dev/null +++ b/src/coreComponents/schema/docs/SourceFluxStatistics.rst @@ -0,0 +1,16 @@ + + +============== ================== ======== ===================================================================================================================================================================================== +Name Type Default Description +============== ================== ======== ===================================================================================================================================================================================== +flowSolverName groupNameRef required Name of the flow solver +fluxNames groupNameRef_array {*} Name(s) array of the SourceFlux(s) for which we want the statistics. Use "*" to target all SourceFlux. +logLevel integer 0 | Log level + | - Log Level 1 outputs the sum of all SourceFlux(s) produced rate & mass, + | - Log Level 2 details values for each SourceFlux, + | - Log Level 3 details values for each region. +name groupName required A name is required for any non-unique nodes +writeCSV integer 0 Write statistics into a CSV file +============== ================== ======== ===================================================================================================================================================================================== + + diff --git a/src/coreComponents/schema/docs/SourceFluxStatistics_other.rst b/src/coreComponents/schema/docs/SourceFluxStatistics_other.rst new file mode 100644 index 00000000000..adf1c1b8aec --- /dev/null +++ b/src/coreComponents/schema/docs/SourceFluxStatistics_other.rst @@ -0,0 +1,9 @@ + + +==== ==== ============================ +Name Type Description +==== ==== ============================ + (no documentation available) +==== ==== ============================ + + diff --git a/src/coreComponents/schema/docs/Tasks.rst b/src/coreComponents/schema/docs/Tasks.rst index 95301497c5c..3a3a3e06a9f 100644 --- a/src/coreComponents/schema/docs/Tasks.rst +++ b/src/coreComponents/schema/docs/Tasks.rst @@ -15,6 +15,7 @@ SinglePhaseReservoirPoromechanicsInitialization node :ref:`X SinglePhaseStatistics node :ref:`XML_SinglePhaseStatistics` SolidMechanicsStateReset node :ref:`XML_SolidMechanicsStateReset` SolidMechanicsStatistics node :ref:`XML_SolidMechanicsStatistics` +SourceFluxStatistics node :ref:`XML_SourceFluxStatistics` TriaxialDriver node :ref:`XML_TriaxialDriver` =========================================================== ==== ======= ====================================================================== diff --git a/src/coreComponents/schema/docs/Tasks_other.rst b/src/coreComponents/schema/docs/Tasks_other.rst index dca37a26a99..4ee34a22e81 100644 --- a/src/coreComponents/schema/docs/Tasks_other.rst +++ b/src/coreComponents/schema/docs/Tasks_other.rst @@ -15,6 +15,7 @@ SinglePhaseReservoirPoromechanicsInitialization node :ref:`DATASTRUC SinglePhaseStatistics node :ref:`DATASTRUCTURE_SinglePhaseStatistics` SolidMechanicsStateReset node :ref:`DATASTRUCTURE_SolidMechanicsStateReset` SolidMechanicsStatistics node :ref:`DATASTRUCTURE_SolidMechanicsStatistics` +SourceFluxStatistics node :ref:`DATASTRUCTURE_SourceFluxStatistics` TriaxialDriver node :ref:`DATASTRUCTURE_TriaxialDriver` =========================================================== ==== ================================================================================ diff --git a/src/coreComponents/schema/schema.xsd b/src/coreComponents/schema/schema.xsd index f7b331fabf4..2180d14bcc5 100644 --- a/src/coreComponents/schema/schema.xsd +++ b/src/coreComponents/schema/schema.xsd @@ -511,6 +511,10 @@ + + + + @@ -1245,7 +1249,7 @@ This keyword is ignored for single-phase flow simulations--> - + @@ -1253,7 +1257,7 @@ This keyword is ignored for single-phase flow simulations--> - + @@ -2424,7 +2428,7 @@ the relative residual norm satisfies: - + @@ -2494,7 +2498,7 @@ the relative residual norm satisfies: - + @@ -3117,7 +3121,11 @@ Local - Add stabilization only to interiors of macro elements.--> - + @@ -3147,7 +3155,11 @@ Local - Add stabilization only to interiors of macro elements.--> - + @@ -3263,7 +3275,11 @@ Local - Add stabilization only to interiors of macro elements.--> - + @@ -3616,6 +3632,7 @@ Local - Add stabilization only to interiors of macro elements.--> + @@ -3781,6 +3798,21 @@ Local - Add stabilization only to interiors of macro elements.--> + + + + + + + + + + + + diff --git a/src/coreComponents/schema/schema.xsd.other b/src/coreComponents/schema/schema.xsd.other index b337dd2d503..e903f4ea2fa 100644 --- a/src/coreComponents/schema/schema.xsd.other +++ b/src/coreComponents/schema/schema.xsd.other @@ -1305,6 +1305,7 @@ + @@ -1320,6 +1321,7 @@ + diff --git a/src/coreComponents/unitTests/CMakeLists.txt b/src/coreComponents/unitTests/CMakeLists.txt index a3d92867343..09542ae3472 100644 --- a/src/coreComponents/unitTests/CMakeLists.txt +++ b/src/coreComponents/unitTests/CMakeLists.txt @@ -10,5 +10,6 @@ add_subdirectory( meshTests ) add_subdirectory( finiteVolumeTests ) add_subdirectory( fileIOTests ) add_subdirectory( fluidFlowTests ) +add_subdirectory( testingUtilities ) add_subdirectory( wellsTests ) -add_subdirectory( wavePropagationTests ) \ No newline at end of file +add_subdirectory( wavePropagationTests ) diff --git a/src/coreComponents/unitTests/fluidFlowTests/CMakeLists.txt b/src/coreComponents/unitTests/fluidFlowTests/CMakeLists.txt index 7abc6ca2de1..2abdc8e16cf 100644 --- a/src/coreComponents/unitTests/fluidFlowTests/CMakeLists.txt +++ b/src/coreComponents/unitTests/fluidFlowTests/CMakeLists.txt @@ -1,10 +1,9 @@ # Specify list of tests - -set(gtest_geosx_tests - testSinglePhaseBaseKernels.cpp - testThermalCompMultiphaseFlow.cpp - testThermalSinglePhaseFlow.cpp - ) +set( gtest_geosx_tests + testSinglePhaseBaseKernels.cpp + testThermalCompMultiphaseFlow.cpp + testThermalSinglePhaseFlow.cpp + testFlowStatistics.cpp ) set( dependencyList ${parallelDeps} gtest ) @@ -14,6 +13,7 @@ else() list( APPEND dependencyList ${geosx_core_libs} ) endif() +set( dependencyList testingUtilities ) if( ENABLE_PVTPackage ) list( APPEND gtest_geosx_tests diff --git a/src/coreComponents/unitTests/fluidFlowTests/testFlowStatistics.cpp b/src/coreComponents/unitTests/fluidFlowTests/testFlowStatistics.cpp new file mode 100644 index 00000000000..12de6b8f655 --- /dev/null +++ b/src/coreComponents/unitTests/fluidFlowTests/testFlowStatistics.cpp @@ -0,0 +1,1146 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2018-2020 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2020 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2020 TotalEnergies + * Copyright (c) 2019- GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +#include "unitTests/fluidFlowTests/testCompFlowUtils.hpp" +#include "unitTests/testingUtilities/TestingTasks.hpp" +#include "mainInterface/initialization.hpp" +#include "mainInterface/GeosxState.hpp" +#include "fieldSpecification/SourceFluxStatistics.hpp" +#include "physicsSolvers/fluidFlow/SinglePhaseStatistics.hpp" + +#include + + +using namespace geos; +using namespace geos::dataRepository; +using namespace geos::testing; + +CommandLineOptions g_commandLineOptions; + + +//////////////////////////////// Test base utilities //////////////////////////////// + + +/** + * @brief this struct is used to provide the input data of each flow tests + */ +struct TestInputs +{ + string xmlInput; + std::map< string, string > tableFiles; + + string sourceFluxName; + string sinkFluxName; + string timeStepCheckerPath; + string timeStepFluxStatsPath; + string wholeSimFluxStatsPath; + string flowSolverPath; + + // rates for each timesteps, for each phases + array2d< real64 > sourceRates; + array2d< real64 > sinkRates; + + // parameters for precomputing results + real64 dt; + real64 sourceRateFactor; + real64 sinkRateFactor; + integer sourceElementsCount; + integer sinkElementsCount; + + /// In order to be sure that sub-timestepping is supported, requires the test to have at least one sub-timestep. + /// At least one simulation should test the timestep cuts ! + integer requiredSubTimeStep = 0; +}; + +/** + * @brief this struct computes from the test inputs the values to expect from the simulation. + */ +struct TestSet +{ + TestInputs const inputs; + + integer timestepCount; + integer totalElementsCount; + integer phaseCount; + + // stats for each timesteps, for each phases + array2d< real64 > sourceRates; + array2d< real64 > sinkRates; + array2d< real64 > sourceMassProd; + array2d< real64 > sinkMassProd; + array2d< real64 > massDeltas; + + // overall simulation stats for each phases + array1d< real64 > sourceMeanRate; + array1d< real64 > sinkMeanRate; + array1d< real64 > totalSourceMassProd; + array1d< real64 > totalSinkMassProd; + array1d< real64 > totalMassProd; + array1d< real64 > totalMeanRate; + + /** + * @brief Compute the expected statistics set for the tested simulation. + * @param inputParams the test simulation input parameters + */ + TestSet( TestInputs const & inputParams ): + inputs( inputParams ) + { + // tables must provide the same timestep & phase rates + EXPECT_EQ( inputs.sourceRates.size( 0 ), inputs.sinkRates.size( 0 )); + EXPECT_EQ( inputs.sourceRates.size( 1 ), inputs.sinkRates.size( 1 )); + + timestepCount = inputs.sourceRates.size( 0 ); + phaseCount = inputs.sourceRates.size( 1 ); + totalElementsCount = inputs.sourceElementsCount + inputs.sinkElementsCount; + + sourceRates.resize( timestepCount, phaseCount ); + sinkRates.resize( timestepCount, phaseCount ); + sourceMassProd.resize( timestepCount, phaseCount ); + sinkMassProd.resize( timestepCount, phaseCount ); + massDeltas.resize( timestepCount, phaseCount ); + sourceMeanRate.resize( phaseCount ); + sinkMeanRate.resize( phaseCount ); + totalSourceMassProd.resize( phaseCount ); + totalSinkMassProd.resize( phaseCount ); + totalMassProd.resize( phaseCount ); + totalMeanRate.resize( phaseCount ); + + for( integer ip = 0; ip < phaseCount; ++ip ) + { + for( integer timestepId = 0; timestepId < timestepCount; ++timestepId ) + { + // mass production / injection calculation + sourceRates[timestepId][ip] = inputs.sourceRates[timestepId][ip] * inputs.sourceRateFactor; + sourceMassProd[timestepId][ip] = inputs.sourceRates[timestepId][ip] * inputs.dt * inputs.sourceRateFactor; + totalSourceMassProd[ip] += sourceMassProd[timestepId][ip]; + sinkRates[timestepId][ip] = inputs.sinkRates[timestepId][ip] * inputs.sinkRateFactor; + sinkMassProd[timestepId][ip] = inputs.sinkRates[timestepId][ip] * inputs.dt * inputs.sinkRateFactor; + massDeltas[timestepId][ip] = -( sourceMassProd[timestepId][ip] + sinkMassProd[timestepId][ip] ); + totalSinkMassProd[ip] += sinkMassProd[timestepId][ip]; + // rates accumulations + sourceMeanRate[ip] += inputs.sourceRates[timestepId][ip] * inputs.sourceRateFactor; + sinkMeanRate[ip] += inputs.sinkRates[timestepId][ip] * inputs.sinkRateFactor; + } + // mean rates calculation + real64 const ratesMeanDivisor = 1.0 / double( timestepCount - 1 ); + sourceMeanRate[ip] *= ratesMeanDivisor; + sinkMeanRate[ip] *= ratesMeanDivisor; + // totals + totalMassProd[ip] = totalSinkMassProd[ip] + totalSourceMassProd[ip]; + totalMeanRate[ip] = sinkMeanRate[ip] + sourceMeanRate[ip]; + } + } +}; + + +class FlowStatisticsTest : public ::testing::Test +{ +public: + + void writeTableFiles( std::map< string, string > const & files ) + { + for( auto const & [fileName, content] : files ) + { + std::ofstream os( fileName ); + ASSERT_TRUE( os.is_open() ); + os << content; + os.close(); + + m_tableFileNames.push_back( fileName ); + } + } + + void TearDown() override + { + // removing temp table files + for( string const & fileName : m_tableFileNames ) + { + ASSERT_TRUE( std::remove( fileName.c_str() ) == 0 ); + } + m_tableFileNames.clear(); + } + +private: + std::vector< string > m_tableFileNames; +}; + + +void setRateTable( array2d< real64 > & rateTable, std::initializer_list< std::initializer_list< real64 > > timestepPhaseValues ) +{ + rateTable.resize( timestepPhaseValues.size(), timestepPhaseValues.begin()->size() ); + integer timestepId = 0; + for( auto const & phaseValues : timestepPhaseValues ) + { + integer ip = 0; + for( auto const & phaseValue : phaseValues ) + { + rateTable[timestepId][ip++] = phaseValue; + } + ++timestepId; + } +} + +real64 getTotalFluidMass( ProblemManager & problem, string_view flowSolverPath ) +{ + real64 totalMass = 0.0; + SolverBase const & solver = problem.getGroupByPath< SolverBase >( string( flowSolverPath ) ); + solver.forDiscretizationOnMeshTargets( problem.getDomainPartition().getMeshBodies(), + [&] ( string const &, + MeshLevel & mesh, + arrayView1d< string const > const & ) + { + mesh.getElemManager().forElementRegions( [&]( ElementRegionBase & region ) + { + SinglePhaseStatistics::RegionStatistics & regionStats = region.getReference< SinglePhaseStatistics::RegionStatistics >( + SinglePhaseStatistics::viewKeyStruct::regionStatisticsString() ); + + totalMass += regionStats.totalMass; + } ); + } ); + return totalMass; +} + + +/** + * @brief Verification that the source flux statistics are correct for the current timestep + * @param expectedMasses the expected mass values per phase + * @param expectedRates the expected rate values per phase + * @param expectedElementCount the number of expected targeted elements + * @param stats the timestep stats + * @param context a context string to provide in any error message. + */ +void checkFluxStats( arraySlice1d< real64 > const & expectedMasses, + arraySlice1d< real64 > const & expectedRates, + integer const expectedElementCount, + SourceFluxStatsAggregator::WrappedStats const & stats, + string_view context ) +{ + for( int ip = 0; ip < stats.stats().getPhaseCount(); ++ip ) + { + EXPECT_DOUBLE_EQ( stats.stats().m_producedMass[ip], expectedMasses[ip] ) << GEOS_FMT( "The flux named '{}' did not produce the expected mass ({}, phase = {}).", + stats.getFluxName(), context, ip ); + EXPECT_DOUBLE_EQ( stats.stats().m_productionRate[ip], expectedRates[ip] ) << GEOS_FMT( "The flux named '{}' did not produce at the expected rate ({}, phase = {}).", + stats.getFluxName(), context, ip ); + } + EXPECT_DOUBLE_EQ( stats.stats().m_elementCount, expectedElementCount ) << GEOS_FMT( "The flux named '{}' did not produce in the expected elements ({}).", + stats.getFluxName(), context ); +} + +/** + * @brief Verification that the source flux statistics are correct over the whole simulation + * @param problem the simulation ProblemManager + * @param testSet the simulation TestSet + */ +void checkWholeSimFluxStats( ProblemManager & problem, TestSet const & testSet ) +{ + DomainPartition & domain = problem.getDomainPartition(); + SourceFluxStatsAggregator & wholeSimStats = + problem.getGroupByPath< SourceFluxStatsAggregator >( testSet.inputs.wholeSimFluxStatsPath ); + + wholeSimStats.forMeshLevelStatsWrapper( domain, + [&] ( MeshLevel & meshLevel, + SourceFluxStatsAggregator::WrappedStats & meshLevelStats ) + { + wholeSimStats.forAllFluxStatsWrappers( meshLevel, + [&] ( MeshLevel &, + SourceFluxStatsAggregator::WrappedStats & fluxStats ) + { + if( fluxStats.getFluxName() == testSet.inputs.sourceFluxName ) + { + checkFluxStats( testSet.totalSourceMassProd, + testSet.sourceMeanRate, + testSet.inputs.sourceElementsCount, + fluxStats, "over whole simulation" ); + } + else if( fluxStats.getFluxName() == testSet.inputs.sinkFluxName ) + { + checkFluxStats( testSet.totalSinkMassProd, + testSet.sinkMeanRate, + testSet.inputs.sinkElementsCount, + fluxStats, "over whole simulation" ); + } + else + { + FAIL() << "Unexpected SourceFlux found!"; + } + } ); + + for( int ip = 0; ip < meshLevelStats.stats().getPhaseCount(); ++ip ) + { + EXPECT_DOUBLE_EQ( meshLevelStats.stats().m_producedMass[ip], testSet.totalMassProd[ip] ) << "The fluxes did not produce the expected total mass (over whole simulation, phase = " << ip << ")."; + EXPECT_DOUBLE_EQ( meshLevelStats.stats().m_productionRate[ip], testSet.totalMeanRate[ip] ) << "The fluxes did not produce at the expected rate (over whole simulation, phase = " << ip << ")."; + } + } ); +} + +/** + * @brief Verification that the source flux statistics are correct for a given timestep + * @param problem the simulation ProblemManager + * @param testSet the simulation TestSet + * @param time_n the current timestep start + * @param timestepId the current timestep id (= cycle) + */ +void checkTimeStepFluxStats( ProblemManager & problem, TestSet const & testSet, + real64 const time_n, integer const timestepId ) +{ + DomainPartition & domain = problem.getDomainPartition(); + SourceFluxStatsAggregator & timestepStats = + problem.getGroupByPath< SourceFluxStatsAggregator >( testSet.inputs.timeStepFluxStatsPath ); + + timestepStats.forMeshLevelStatsWrapper( domain, + [&] ( MeshLevel & meshLevel, + SourceFluxStatsAggregator::WrappedStats & ) + { + timestepStats.forAllFluxStatsWrappers( meshLevel, + [&] ( MeshLevel &, + SourceFluxStatsAggregator::WrappedStats & fluxStats ) + { + if( fluxStats.getFluxName() == testSet.inputs.sourceFluxName ) + { + checkFluxStats( testSet.sourceMassProd[timestepId], + testSet.sourceRates[timestepId], + testSet.inputs.sourceElementsCount, + fluxStats, GEOS_FMT( "for timestep at t = {} s", time_n ) ); + } + else if( fluxStats.getFluxName() == testSet.inputs.sinkFluxName ) + { + checkFluxStats( testSet.sinkMassProd[timestepId], + testSet.sinkRates[timestepId], + testSet.inputs.sinkElementsCount, + fluxStats, GEOS_FMT( "for timestep at t = {} s", time_n ) ); + } + else + { + FAIL() << "Unexpected SourceFlux found!"; + } + } ); + } ); +} + +void checkTimeStepStats( TestSet const & testSet, + real64 const time_n, + integer const timestepId ) +{ + EXPECT_LT( timestepId, testSet.timestepCount ) << GEOS_FMT( "The tested time-step count were higher than expected (t = {} s).", + time_n ); +} + +void checkWholeSimTimeStepStats( ProblemManager & problem, + TestSet const & testSet, + TimeStepChecker const & timeStepChecker ) +{ + EXPECT_EQ( timeStepChecker.getTestedTimeStepCount(), testSet.timestepCount ) << "The tested time-step were different than expected."; + + SolverBase const & solver = problem.getGroupByPath< SolverBase >( testSet.inputs.flowSolverPath ); + SolverStatistics const & solverStats = solver.getSolverStatistics(); + EXPECT_GE( solverStats.getNumTimeStepCuts(), testSet.inputs.requiredSubTimeStep ) << "The test did not encountered any timestep cut, but were expected to. " + "Consider adapting the simulation so a timestep cut occurs to check they work as expected."; +} + + +//////////////////////////////// SinglePhase Flux Statistics Test //////////////////////////////// +namespace SinglePhaseFluxStatisticsTest +{ + + +TestSet getTestSet() +{ + TestInputs testInputs; + + testInputs.xmlInput = + R"xml( + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +)xml"; + + + testInputs.sourceFluxName = "sourceFlux"; + testInputs.sinkFluxName = "sinkFlux"; + testInputs.timeStepCheckerPath = "/Tasks/timeStepChecker"; + testInputs.timeStepFluxStatsPath = "/Tasks/timeStepFluxStats"; + testInputs.wholeSimFluxStatsPath = "/Tasks/wholeSimFluxStats"; + testInputs.flowSolverPath = "/Solvers/testSolver"; + + testInputs.dt = 500.0; + testInputs.sourceElementsCount = 2; + testInputs.sinkElementsCount = 5; + + // FluxRate table from 0.0s to 5000.0s + setRateTable( testInputs.sourceRates, + { { 0.000 }, + { 0.000 }, + { 0.767 }, + { 0.894 }, + { 0.561 }, + { 0.234 }, + { 0.194 }, + { 0.178 }, + { 0.162 }, + { 0.059 }, + { 0.000 } } ); + testInputs.sinkRates=testInputs.sourceRates; + + // sink is 3x source production + testInputs.sourceRateFactor = -1.0; + testInputs.sinkRateFactor = 3.0; + + return TestSet( testInputs ); +} + +TEST_F( FlowStatisticsTest, checkSinglePhaseFluxStatistics ) +{ + TestSet const testSet = getTestSet(); + + GeosxState state( std::make_unique< CommandLineOptions >( g_commandLineOptions ) ); + ProblemManager & problem = state.getProblemManager(); + + setupProblemFromXML( problem, testSet.inputs.xmlInput.data() ); + + real64 firstMass; + + TimeStepChecker & timeStepChecker = problem.getGroupByPath< TimeStepChecker >( testSet.inputs.timeStepCheckerPath ); + timeStepChecker.setTimeStepCheckingFunction( [&]( real64 const time_n ) + { + integer const timestepId = timeStepChecker.getTestedTimeStepCount(); + checkTimeStepStats( testSet, time_n, timestepId ); + checkTimeStepFluxStats( problem, testSet, time_n, timestepId ); + + static bool passedFirstTimeStep = false; + if( !passedFirstTimeStep ) + { + passedFirstTimeStep = true; + firstMass = getTotalFluidMass( problem, testSet.inputs.flowSolverPath ); + } + } ); + + // run simulation + EXPECT_FALSE( problem.runSimulation() ) << "Simulation exited early."; + + checkWholeSimFluxStats( problem, testSet ); + checkWholeSimTimeStepStats( problem, testSet, timeStepChecker ); + + // check singlephasestatistics results + real64 const lastMass = getTotalFluidMass( problem, testSet.inputs.flowSolverPath ); + real64 const massDiffTol = 1e-7; + EXPECT_NEAR( lastMass - firstMass, + -testSet.totalMassProd[0], + massDiffTol * std::abs( testSet.totalMassProd[0] ) ) << GEOS_FMT( "{} total mass difference from start to end is not consistent with fluxes production.", + SinglePhaseStatistics::catalogName() ); +} + + +} /* namespace SinglePhaseFluxStatisticsTest */ + + +//////////////////////////////// Multiphase Flux Statistics Test //////////////////////////////// +namespace MultiPhaseFluxStatisticsTestMass +{ + + +TestSet getTestSet() +{ + TestInputs testInputs; + + testInputs.xmlInput = + R"xml( + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +)xml"; + + testInputs.tableFiles["pvtgas.txt"] = "DensityFun SpanWagnerCO2Density 1.5e7 2.5e7 1e5 370.15 400.15 2\n" + "ViscosityFun FenghourCO2Viscosity 1.5e7 2.5e7 1e5 370.15 400.15 2\n"; + + testInputs.tableFiles["pvtliquid.txt"] = "DensityFun EzrokhiBrineDensity 0.1033 -2.2991e-5 -2.3658e-6\n" + "ViscosityFun EzrokhiBrineViscosity 0 0 0\n"; + + testInputs.tableFiles["co2flash.txt"] = "FlashModel CO2Solubility 1.5e7 2.5e7 1e5 370.15 400.15 2 0\n"; + + + testInputs.sourceFluxName = "sourceFlux"; + testInputs.sinkFluxName = "sinkFlux"; + testInputs.timeStepCheckerPath = "/Tasks/timeStepChecker"; + testInputs.timeStepFluxStatsPath = "/Tasks/timeStepFluxStats"; + testInputs.wholeSimFluxStatsPath = "/Tasks/wholeSimFluxStats"; + testInputs.flowSolverPath = "/Solvers/testSolver"; + + testInputs.dt = 500.0; + testInputs.sourceElementsCount = 1; + testInputs.sinkElementsCount = 1; + + // FluxInjectionRate & FluxProductionRate table from 0.0s to 5000.0s + setRateTable( testInputs.sourceRates, + { { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.267, 0.0 }, + { 0.561, 0.0 }, + { 0.194, 0.0 }, + { 0.102, 0.0 }, + { 0.059, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 } } ); + setRateTable( testInputs.sinkRates, + { { 0.0, 0.000 }, + { 0.0, 0.000 }, + { 0.0, 0.003 }, + { 0.0, 0.062 }, + { 0.0, 0.121 }, + { 0.0, 0.427 }, + { 0.0, 0.502 }, + { 0.0, 0.199 }, + { 0.0, 0.083 }, + { 0.0, 0.027 }, + { 0.0, 0.000 } } ); + + testInputs.sourceRateFactor = -44e-3; + testInputs.sinkRateFactor = 18e-3; + + return TestSet( testInputs ); +} + + +TEST_F( FlowStatisticsTest, checkMultiPhaseFluxStatisticsMass ) +{ + TestSet const testSet = getTestSet(); + writeTableFiles( testSet.inputs.tableFiles ); + + GeosxState state( std::make_unique< CommandLineOptions >( g_commandLineOptions ) ); + ProblemManager & problem = state.getProblemManager(); + + setupProblemFromXML( problem, testSet.inputs.xmlInput.data() ); + + TimeStepChecker & timeStepChecker = problem.getGroupByPath< TimeStepChecker >( testSet.inputs.timeStepCheckerPath ); + timeStepChecker.setTimeStepCheckingFunction( [&]( real64 const time_n ) + { + integer const timestepId = timeStepChecker.getTestedTimeStepCount(); + checkTimeStepStats( testSet, time_n, timestepId ); + checkTimeStepFluxStats( problem, testSet, time_n, timestepId ); + } ); + + // run simulation + EXPECT_FALSE( problem.runSimulation() ) << "Simulation exited early."; + + checkWholeSimFluxStats( problem, testSet ); + checkWholeSimTimeStepStats( problem, testSet, timeStepChecker ); +} + + +} /* namespace MultiPhaseFluxStatisticsTest */ + + +//////////////////////////////// Multiphase Flux Statistics Test //////////////////////////////// +namespace MultiPhaseFluxStatisticsTestMol +{ + + +TestSet getTestSet() +{ + TestInputs testInputs; + + testInputs.xmlInput = + R"xml( + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +)xml"; + + testInputs.tableFiles["pvtgas.txt"] = "DensityFun SpanWagnerCO2Density 1.5e7 2.5e7 1e5 370.15 400.15 2\n" + "ViscosityFun FenghourCO2Viscosity 1.5e7 2.5e7 1e5 370.15 400.15 2\n"; + + testInputs.tableFiles["pvtliquid.txt"] = "DensityFun EzrokhiBrineDensity 0.1033 -2.2991e-5 -2.3658e-6\n" + "ViscosityFun EzrokhiBrineViscosity 0 0 0\n"; + + testInputs.tableFiles["co2flash.txt"] = "FlashModel CO2Solubility 1.5e7 2.5e7 1e5 370.15 400.15 2 0\n"; + + + testInputs.sourceFluxName = "sourceFlux"; + testInputs.sinkFluxName = "sinkFlux"; + testInputs.timeStepCheckerPath = "/Tasks/timeStepChecker"; + testInputs.timeStepFluxStatsPath = "/Tasks/timeStepFluxStats"; + testInputs.wholeSimFluxStatsPath = "/Tasks/wholeSimFluxStats"; + testInputs.flowSolverPath = "/Solvers/testSolver"; + + testInputs.dt = 500.0; + testInputs.sourceElementsCount = 1; + testInputs.sinkElementsCount = 1; + + // FluxInjectionRate & FluxProductionRate table from 0.0s to 5000.0s + setRateTable( testInputs.sourceRates, + { { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.267, 0.0 }, + { 0.561, 0.0 }, + { 0.194, 0.0 }, + { 0.102, 0.0 }, + { 0.059, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 }, + { 0.000, 0.0 } } ); + setRateTable( testInputs.sinkRates, + { { 0.0, 0.000 }, + { 0.0, 0.000 }, + { 0.0, 0.003 }, + { 0.0, 0.062 }, + { 0.0, 0.121 }, + { 0.0, 0.427 }, + { 0.0, 0.502 }, + { 0.0, 0.199 }, + { 0.0, 0.083 }, + { 0.0, 0.027 }, + { 0.0, 0.000 } } ); + + // scale is set to high values to make the solver generate timestep cuts + testInputs.sourceRateFactor = -8.0; + testInputs.sinkRateFactor = 8.0; + + // this simulation is set-up to have at least one timestep cut. + testInputs.requiredSubTimeStep = 2; + + return TestSet( testInputs ); +} + + +TEST_F( FlowStatisticsTest, checkMultiPhaseFluxStatisticsMol ) +{ + TestSet const testSet = getTestSet(); + writeTableFiles( testSet.inputs.tableFiles ); + + GeosxState state( std::make_unique< CommandLineOptions >( g_commandLineOptions ) ); + ProblemManager & problem = state.getProblemManager(); + + setupProblemFromXML( problem, testSet.inputs.xmlInput.data() ); + + TimeStepChecker & timeStepChecker = problem.getGroupByPath< TimeStepChecker >( testSet.inputs.timeStepCheckerPath ); + timeStepChecker.setTimeStepCheckingFunction( [&]( real64 const time_n ) + { + integer const timestepId = timeStepChecker.getTestedTimeStepCount(); + checkTimeStepStats( testSet, time_n, timestepId ); + checkTimeStepFluxStats( problem, testSet, time_n, timestepId ); + } ); + + // run simulation + EXPECT_FALSE( problem.runSimulation() ) << "Simulation exited early."; + + checkWholeSimFluxStats( problem, testSet ); + checkWholeSimTimeStepStats( problem, testSet, timeStepChecker ); +} + + +} /* namespace MultiPhaseFluxStatisticsTest */ + + +//////////////////////////////// Main //////////////////////////////// + + +int main( int argc, char * * argv ) +{ + ::testing::InitGoogleTest( &argc, argv ); + g_commandLineOptions = *geos::basicSetup( argc, argv ); + int const result = RUN_ALL_TESTS(); + geos::basicCleanup(); + return result; +} diff --git a/src/coreComponents/unitTests/testingUtilities/CMakeLists.txt b/src/coreComponents/unitTests/testingUtilities/CMakeLists.txt new file mode 100644 index 00000000000..f6d8d8632b0 --- /dev/null +++ b/src/coreComponents/unitTests/testingUtilities/CMakeLists.txt @@ -0,0 +1,33 @@ +# +# Specify all headers +# +set( testingUtilities_headers + TestingTasks.hpp + ) + +# +# Specify all sources +# +set( testingUtilities_sources + TestingTasks.cpp + ) + +# +# Specify all dependencies +# +set( dependencyList ${parallelDeps} gtest ) + +if ( GEOSX_BUILD_SHARED_LIBS ) + list( APPEND dependencyList geosx_core ) +else() + list( APPEND dependencyList ${geosx_core_libs} ) +endif() + +blt_add_library( NAME testingUtilities + SOURCES ${testingUtilities_sources} + HEADERS ${testingUtilities_headers} + DEPENDS_ON ${dependencyList} + OBJECT ${GEOSX_BUILD_OBJ_LIBS} + ) + +target_include_directories( testingUtilities PUBLIC ${CMAKE_SOURCE_DIR}/coreComponents ) diff --git a/src/coreComponents/unitTests/testingUtilities/TestingTasks.cpp b/src/coreComponents/unitTests/testingUtilities/TestingTasks.cpp new file mode 100644 index 00000000000..b948a406071 --- /dev/null +++ b/src/coreComponents/unitTests/testingUtilities/TestingTasks.cpp @@ -0,0 +1,50 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2018-2020 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2020 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2020 TotalEnergies + * Copyright (c) 2019- GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file TestingTasks.cpp + */ + +#include "TestingTasks.hpp" + +namespace geos +{ +namespace testing +{ + + +TimeStepChecker::TimeStepChecker( string const & name, Group * const parent ): + TaskBase( name, parent ) +{} + +bool TimeStepChecker::execute( real64 const time_n, + real64 const GEOS_UNUSED_PARAM( dt ), + integer const GEOS_UNUSED_PARAM( cycleNumber ), + integer const GEOS_UNUSED_PARAM( eventCounter ), + real64 const GEOS_UNUSED_PARAM( eventProgress ), + DomainPartition & GEOS_UNUSED_PARAM( domain ) ) +{ + EXPECT_TRUE( m_checkTimeStepFunction ); + m_checkTimeStepFunction( time_n ); + + ++m_timestepId; + return false; +} + +REGISTER_CATALOG_ENTRY( TaskBase, TimeStepChecker, string const &, geos::dataRepository::Group * const ) + + +} // namespace testing + +} // namespace geos diff --git a/src/coreComponents/unitTests/testingUtilities/TestingTasks.hpp b/src/coreComponents/unitTests/testingUtilities/TestingTasks.hpp new file mode 100644 index 00000000000..f52a9e54834 --- /dev/null +++ b/src/coreComponents/unitTests/testingUtilities/TestingTasks.hpp @@ -0,0 +1,80 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2018-2020 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2020 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2020 TotalEnergies + * Copyright (c) 2019- GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/* + * @file TestingTasks.hpp + */ + +#ifndef GEOS_EVENTS_TASKS_TESTINGTASKS_HPP_ +#define GEOS_EVENTS_TASKS_TESTINGTASKS_HPP_ + +#include "events/tasks/TaskBase.hpp" +#include "codingUtilities/UnitTestUtilities.hpp" + +namespace geos +{ +namespace testing +{ + + +/** + * @brief This Task allows to do checks each timestep during the simulation by calling a provided functor. + * To be executed, it must be - as every task - declared in the provided xmlInput within the Tasks node and called by an even. + * As this Group is designed for developpers and not user oriented, it should not appear in the documentation nor in the xsd. + * Question: could this task be used in the integratedTests ? + */ +class TimeStepChecker : public TaskBase +{ +public: + /** + * @brief Construct a new TimeStepChecker Task + * @param name name in xsd + * @param parent parent group in hierarchy + */ + TimeStepChecker( string const & name, Group * const parent ); + + /** + * @brief Set the functor that must be called each time this Task is executed. + * @param func the functor to execute + */ + void setTimeStepCheckingFunction( std::function< void(real64) > func ) + { m_checkTimeStepFunction = std::function< void(real64) >( func ); } + + /** + * @return Get the tested time-step count + */ + integer getTestedTimeStepCount() const + { return m_timestepId; } + + /** + * @brief Catalog name interface + * @return This type's catalog name + */ + static string catalogName() { return "TimeStepChecker"; } + + virtual bool execute( real64 time_n, real64 dt, integer cycleNumber, + integer eventCounter, real64 eventProgress, + DomainPartition & domain ); + +private: + std::function< void(real64) > m_checkTimeStepFunction; + int m_timestepId = 0; +}; + + +} // namespace testing + +} // namespace geos + +#endif //GEOS_EVENTS_TASKS_TESTINGTASKS_HPP_ diff --git a/src/docs/sphinx/CompleteXMLSchema.rst b/src/docs/sphinx/CompleteXMLSchema.rst index 7ab5a264d15..46a0ade984c 100644 --- a/src/docs/sphinx/CompleteXMLSchema.rst +++ b/src/docs/sphinx/CompleteXMLSchema.rst @@ -1235,6 +1235,13 @@ Element: SourceFlux .. include:: ../../coreComponents/schema/docs/SourceFlux.rst +.. _XML_SourceFluxStatistics: + +Element: SourceFluxStatistics +============================= +.. include:: ../../coreComponents/schema/docs/SourceFluxStatistics.rst + + .. _XML_SurfaceElementRegion: Element: SurfaceElementRegion @@ -2695,6 +2702,13 @@ Datastructure: SourceFlux .. include:: ../../coreComponents/schema/docs/SourceFlux_other.rst +.. _DATASTRUCTURE_SourceFluxStatistics: + +Datastructure: SourceFluxStatistics +=================================== +.. include:: ../../coreComponents/schema/docs/SourceFluxStatistics_other.rst + + .. _DATASTRUCTURE_SurfaceElementRegion: Datastructure: SurfaceElementRegion