diff --git a/inputFiles/singlePhaseFlow/sourceFlux_2d.xml b/inputFiles/singlePhaseFlow/sourceFlux_2d.xml index 4a0d5cf3a28..9452ee415ff 100644 --- a/inputFiles/singlePhaseFlow/sourceFlux_2d.xml +++ b/inputFiles/singlePhaseFlow/sourceFlux_2d.xml @@ -145,31 +145,12 @@ scale="5e6"/> - - - - + scale="{ 2.0e-22, 2.0e-22, 2.0e-22 }"/> #include "dataRepository/xmlWrapper.hpp" +#include "codingUtilities/RTTypes.hpp" #include "common/format/EnumStrings.hpp" using namespace geos; @@ -443,6 +444,10 @@ TEST( testXmlWrapper, testGroupNamesArrayFormats ) GroupNameTest( groupNameRefArrayRegex, "{path/to/resource*}" ), GroupNameTest( groupNameRefArrayRegex, "{[arrayElement]}" ), GroupNameTest( groupNameRefArrayRegex, "{name1,name2,name3}" ), + GroupNameTest( groupNameRefArrayRegex, "groupName" ), + GroupNameTest( groupNameRefArrayRegex, " groupName" ), + GroupNameTest( groupNameRefArrayRegex, "name.with-special_chars" ), + GroupNameTest( groupNameRefArrayRegex, "path/to/resource*" ), }; for( GroupNameTest const & input : workingInputs ) { @@ -475,7 +480,6 @@ TEST( testXmlWrapper, testGroupNamesArrayFormats ) GroupNameTest( groupNameRefArrayRegex, "{test\tname}" ), GroupNameTest( groupNameRefArrayRegex, "{test\nname}" ), GroupNameTest( groupNameRefArrayRegex, "{test\rname}" ), - GroupNameTest( groupNameRefArrayRegex, "groupName" ), GroupNameTest( groupNameRefArrayRegex, "{groupName" ), GroupNameTest( groupNameRefArrayRegex, "groupName}" ), GroupNameTest( groupNameRefArrayRegex, "{groupName}} " ), @@ -491,6 +495,9 @@ TEST( testXmlWrapper, testGroupNamesArrayFormats ) GroupNameTest( groupNameRefArrayRegex, "{valuewith,,commas }" ), GroupNameTest( groupNameRefArrayRegex, "{ value with , commas }" ), GroupNameTest( groupNameRefArrayRegex, "{{groupname}}" ), + GroupNameTest( groupNameRefArrayRegex, "group name" ), + GroupNameTest( groupNameRefArrayRegex, "group,name" ), + GroupNameTest( groupNameRefArrayRegex, "group;name" ), }; for( GroupNameTest const & input : erroneousInputs ) { @@ -501,6 +508,53 @@ TEST( testXmlWrapper, testGroupNamesArrayFormats ) } } +TEST( testXmlWrapper, testUnbracketedSingleValueFor1dArrays ) +{ + { + array1d< real64 > array; + EXPECT_NO_THROW( xmlWrapper::stringToInputVariable( array, "1.5", + rtTypes::getTypeRegex< array1d< real64 > >() ) ); + ASSERT_EQ( array.size(), 1 ); + EXPECT_DOUBLE_EQ( array[0], 1.5 ); + } + { + array1d< real64 > array; + EXPECT_NO_THROW( xmlWrapper::stringToInputVariable( array, " -2.3e-4 ", + rtTypes::getTypeRegex< array1d< real64 > >() ) ); + ASSERT_EQ( array.size(), 1 ); + EXPECT_DOUBLE_EQ( array[0], -2.3e-4 ); + } + { + array1d< integer > array; + EXPECT_NO_THROW( xmlWrapper::stringToInputVariable( array, "42", + rtTypes::getTypeRegex< array1d< integer > >() ) ); + ASSERT_EQ( array.size(), 1 ); + EXPECT_EQ( array[0], 42 ); + } + { + array1d< real64 > array; + EXPECT_NO_THROW( xmlWrapper::stringToInputVariable( array, "{ 1, 2, 3 }", + rtTypes::getTypeRegex< array1d< real64 > >() ) ); + ASSERT_EQ( array.size(), 3 ); + EXPECT_DOUBLE_EQ( array[0], 1.0 ); + EXPECT_DOUBLE_EQ( array[1], 2.0 ); + EXPECT_DOUBLE_EQ( array[2], 3.0 ); + } + { + stdVector< std::string > array; + Regex const & arrayRegex = rtTypes::getTypeRegex< string >( rtTypes::CustomTypes::groupNameRefArray ); + EXPECT_NO_THROW( xmlWrapper::stringToInputVariable( array, "myRegion", + arrayRegex ) ); + ASSERT_EQ( array.size(), 1 ); + EXPECT_EQ( array[0], "myRegion" ); + } + { + array2d< real64 > array; + EXPECT_ANY_THROW( xmlWrapper::stringToInputVariable( array, "1.5", + rtTypes::getTypeRegex< array1d< real64 > >() ) ); + } +} + int main( int argc, char * argv[] ) { diff --git a/src/coreComponents/dataRepository/xmlWrapper.cpp b/src/coreComponents/dataRepository/xmlWrapper.cpp index 0c5a469a570..5b172efed0d 100644 --- a/src/coreComponents/dataRepository/xmlWrapper.cpp +++ b/src/coreComponents/dataRepository/xmlWrapper.cpp @@ -136,7 +136,12 @@ void stringToInputVariable( stdVector< std::string > & array, string const & val { validateString( value, regex ); array1d< std::string > tmp; - LvArray::input::stringToArray( tmp, string( stringutilities::trimSpaces( value ) ) ); + string valueTrimmed{ stringutilities::trimSpaces( value ) }; + if( !valueTrimmed.empty() && valueTrimmed.front() != '{' ) + { + valueTrimmed = "{ " + valueTrimmed + " }"; + } + LvArray::input::stringToArray( tmp, valueTrimmed ); array.resize( tmp.size() ); for( localIndex i = 0; i < tmp.size(); ++i ) diff --git a/src/coreComponents/dataRepository/xmlWrapper.hpp b/src/coreComponents/dataRepository/xmlWrapper.hpp index 494eb2db119..081921110f8 100644 --- a/src/coreComponents/dataRepository/xmlWrapper.hpp +++ b/src/coreComponents/dataRepository/xmlWrapper.hpp @@ -391,7 +391,12 @@ stringToInputVariable( Array< T, NDIM, PERMUTATION > & array, string const & val { validateString( value, regex ); - LvArray::input::stringToArray( array, string( stringutilities::trimSpaces( value ) ) ); + string valueTrimmed{ stringutilities::trimSpaces( value ) }; + if( NDIM == 1 && !valueTrimmed.empty() && valueTrimmed.front() != '{' ) + { + valueTrimmed = "{ " + valueTrimmed + " }"; + } + LvArray::input::stringToArray( array, valueTrimmed ); } /** diff --git a/src/coreComponents/fieldSpecification/AquiferBoundaryCondition.cpp b/src/coreComponents/fieldSpecification/AquiferBoundaryCondition.cpp index e420fc303dc..6627fc299e8 100644 --- a/src/coreComponents/fieldSpecification/AquiferBoundaryCondition.cpp +++ b/src/coreComponents/fieldSpecification/AquiferBoundaryCondition.cpp @@ -116,6 +116,8 @@ AquiferBoundaryCondition::AquiferBoundaryCondition( string const & name, Group * void AquiferBoundaryCondition::postInputInitialization() { + FieldSpecification::postInputInitialization(); + GEOS_THROW_IF_LE_MSG( m_permeability, 0.0, "The aquifer permeability cannot be equal to zero or negative", InputError, getDataContext() ); diff --git a/src/coreComponents/fieldSpecification/EquilibriumInitialCondition.cpp b/src/coreComponents/fieldSpecification/EquilibriumInitialCondition.cpp index 8654e47b1b5..d92fe6a1893 100644 --- a/src/coreComponents/fieldSpecification/EquilibriumInitialCondition.cpp +++ b/src/coreComponents/fieldSpecification/EquilibriumInitialCondition.cpp @@ -97,6 +97,7 @@ EquilibriumInitialCondition::EquilibriumInitialCondition( string const & name, G void EquilibriumInitialCondition::postInputInitialization() { + FieldSpecification::postInputInitialization(); FunctionManager const & functionManager = FunctionManager::getInstance(); diff --git a/src/coreComponents/fieldSpecification/FieldSpecification.cpp b/src/coreComponents/fieldSpecification/FieldSpecification.cpp index 74f910d6cb4..7f5c33560d2 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecification.cpp +++ b/src/coreComponents/fieldSpecification/FieldSpecification.cpp @@ -55,19 +55,15 @@ FieldSpecification::FieldSpecification( string const & name, Group * parent ): setDescription( "Direction to apply boundary condition to." ); registerWrapper( viewKeyStruct::functionNameString(), &m_functionName ). - setRTTypeName( rtTypes::CustomTypes::groupNameRef ). - setInputFlag( InputFlags::OPTIONAL ). - setDescription( "Name of function that specifies variation of the boundary condition." ); - - registerWrapper( viewKeyStruct::bcApplicationTableNameString(), &m_bcApplicationFunctionName ). - setRTTypeName( rtTypes::CustomTypes::groupNameRef ). + setRTTypeName( rtTypes::CustomTypes::groupNameRefArray ). setInputFlag( InputFlags::OPTIONAL ). - setDescription( "Name of table that specifies the on/off application of the boundary condition." ); + setDescription( "Name(s) of function(s) that specifies variation of the boundary condition." ); registerWrapper( viewKeyStruct::scaleString(), &m_scale ). setApplyDefaultValue( 0.0 ). setInputFlag( InputFlags::OPTIONAL ). - setDescription( "Apply a scaling factor for the value of the boundary condition." ); + setSizedFromParent( 0 ). + setDescription( "Apply scaling factor(s) for the value(s) of the boundary condition." ); registerWrapper( viewKeyStruct::initialConditionString(), &m_initialCondition ). setApplyDefaultValue( 0 ). @@ -108,6 +104,29 @@ FieldSpecification::getCatalog() } +void FieldSpecification::postInputInitialization() +{ + GEOS_THROW_IF( !m_functionName.empty() && + m_functionName.size() != 1 && + m_functionName.size() != static_cast< string_array::size_type >( m_scale.size() ), + GEOS_FMT ( "Size mismatch: '{}' has {} entries but '{}' has {}. " + "'{}' either must be empty, have a single entry, or be sized exactly like '{}'", + viewKeyStruct::functionNameString(), m_functionName.size(), + viewKeyStruct::scaleString(), m_scale.size(), + viewKeyStruct::functionNameString(), viewKeyStruct::scaleString() ), + InputError, + getDataContext() ); + + if( usesNonScalarValues() ) + { + GEOS_THROW_IF( m_component != -1, + GEOS_FMT ( "'{}' must not be set when '{}' has more than one value.", + viewKeyStruct::componentString(), + viewKeyStruct::scaleString() ), + InputError, + getDataContext() ); + } +} void FieldSpecification::setMeshObjectPath( Group const & meshBodies ) { diff --git a/src/coreComponents/fieldSpecification/FieldSpecification.hpp b/src/coreComponents/fieldSpecification/FieldSpecification.hpp index da0caa4195c..3e52826a26f 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecification.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecification.hpp @@ -129,8 +129,6 @@ class FieldSpecification : public dataRepository::Group constexpr static char const * componentString() { return "component"; } /// @return The key for direction constexpr static char const * directionString() { return "direction"; } - /// @return The key for bcApplicationTableName - constexpr static char const * bcApplicationTableNameString() { return "bcApplicationTableName"; } /// @return The key for scale constexpr static char const * scaleString() { return "scale"; } /// @return The key for functionName @@ -147,10 +145,26 @@ class FieldSpecification : public dataRepository::Group /** * Accessor - * @return const reference to m_function + * @return first entry of m_functionName, or an empty string if empty + * + * @note Legacy scalar accessor. + * Use getFunctionNames() to access the full list of function names when using non-scalar + * field specifications (eg. functionName="{ f1, f2, f3 }") */ string const & getFunctionName() const - { return m_functionName; } + { + static string const emptyName; + return m_functionName.empty() ? emptyName : m_functionName.front(); + } + + /** + * Accessor + * @return const reference to m_functionNames + */ + string_array const & getFunctionNames() const + { + return m_functionName; + } /** * Accessor @@ -210,10 +224,21 @@ class FieldSpecification : public dataRepository::Group /** * Accessor - * @return const m_scale + * @return first entry of m_scale, or 0 if m_scale is empty + * + * @note Legacy scalar accessor. + * Use getScales() to access the full list of scales when using non-scalar + * field specifications (eg. scales="{ 1, 2, 3 }") */ real64 getScale() const - { return m_scale; } + { return m_scale.empty() ? 0.0 : m_scale.front(); } + + /** + * Accessor + * @return const m_scales + */ + arrayView1d< real64 const > getScales() const + { return m_scale.toViewConst(); } /** * Mutator @@ -234,7 +259,27 @@ class FieldSpecification : public dataRepository::Group * @param[in] scale Scaling factor */ void setScale( real64 const & scale ) - { m_scale = scale; } + { + m_scale.resize( 1 ); + m_scale[ 0 ] = scale; + } + + /** + * Mutator + * @brief Set the per-component scale factors + * @param[in] scales The tensor-valued scale + */ + void setScales( array1d< real64 > const & scales ) + { m_scale = scales; } + + /** + * Mutator + * @brief Set the per-component function names + * @param[in] functionNames The per-component function names. Must either be empty, + * have a single entry, or be sized exactly as @p m_scale + */ + void setFunctionNames( string_array const & functionNames ) + { m_functionName = functionNames; } /** * Mutator @@ -265,9 +310,19 @@ class FieldSpecification : public dataRepository::Group MeshObjectPath const & getMeshObjectPaths() const { return *(m_meshObjectPaths.get()); } + /** + * @brief Query whether this field specification uses non-scalar values + * @return True if the field specification provides more than one value + * in the 'scale' or 'functionName' attributes + */ + bool usesNonScalarValues() const + { return m_scale.size() > 1 || m_functionName.size() > 1; } + protected: + virtual void postInputInitialization() override; + private: @@ -294,11 +349,11 @@ class FieldSpecification : public dataRepository::Group /// Whether or not the boundary condition is an initial condition. int m_initialCondition; - /// The name of the function used to generate values for application. - string m_functionName; + /// Name(s) of the function used to generate values for application. + string_array m_functionName; - /// The scale factor to use on the value of the boundary condition. - real64 m_scale; + /// Scale factor(s) to use on the value of the boundary condition. + array1d< real64 > m_scale; /// Time after which the bc is allowed to be applied real64 m_beginTime; @@ -306,9 +361,6 @@ class FieldSpecification : public dataRepository::Group /// Time after which the bc will no longer be applied. real64 m_endTime; - /// The name of a function used to turn on and off the boundary condition. - string m_bcApplicationFunctionName; - /// Enum containing the possible output modes when an error occur SetErrorMode m_emptySetErrorMode; }; diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationImpl.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationImpl.hpp index e23419628f1..8031d121c4d 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationImpl.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationImpl.hpp @@ -314,6 +314,57 @@ class FieldSpecificationImpl arrayView1d< real64 > const & rhsContribution, LAMBDA && lambda ); + /** + * @brief Compute the contributions that will be added/enforced to the right-hand side, + * and collect the corresponding dof numbers + * @tparam FIELD_OP A wrapper struct to define how the boundary condition operates on the variables. + * Either \ref OpEqual or \ref OpAdd. + * @tparam POLICY Execution policy to use when iterating over target set. + * @tparam LAMBDA The type of lambda function passed into the parameter list. + * @param[in] fs The field specification data object. + * @param[in] component The field component to apply the boundary condition to. + * @param[in] scale The scale factor to apply to the boundary condition value. + * @param[in] functionName The name of the function used to evaluate the boundary condition value. + * @param[in] targetSet The set of indices which the boundary condition will be applied. + * @param[in] time The time at which any time dependent functions are to be evaluated as part of the + * application of the boundary condition. + * @param[in] dt time step size which is applied as a factor to bc values + * @param[in] dataGroup The Group that contains the field to apply the boundary condition to. + * @param[in] dofMap The map from the local index of the primary field to the global degree of + * freedom number. + * @param[in] dofRankOffset Offset of dof indices on current rank. + * @param[inout] matrix Local part of the system matrix. + * @param[inout] dof array storing the degrees of freedom of the rhsContribution, to know where + * in the rhs they will be added/enforced + * @param[inout] rhsContribution array storing the values that will be added/enforced to the right-hand side + * @param[in] lambda A lambda function which defines how the value that is passed into the functions + * provided by the FIELD_OP templated type. + * + * This overload behaves like the legacy computeRhsContribution, but takes the component, scale + * and functionName as explicit arguments rather than reading them from the corresponding members. + * This is the variant called when using non-scalar valued field specifications/boundary conditions, + * where each component is applied in turn. + * + * Note that this function only computes the rhs contributions, but does not apply them to the right-hand side. + * The application of these rhs contributions is done in applyBoundaryConditionToSystem. + */ + template< typename FIELD_OP, typename POLICY, typename LAMBDA > + static void + computeRhsContribution( FieldSpecification const & fs, + integer const component, + real64 const scale, + string const & functionName, + SortedArrayView< localIndex const > const & targetSet, + real64 const time, + real64 const dt, + dataRepository::Group const & dataGroup, + arrayView1d< globalIndex const > const & dofMap, + globalIndex const dofRankOffset, + CRSMatrixView< real64, globalIndex const > const & matrix, + arrayView1d< globalIndex > const & dof, + arrayView1d< real64 > const & rhsContribution, + LAMBDA && lambda ); + /** * @brief Function to zero matrix rows to apply boundary conditions * @tparam POLICY the execution policy to use when zeroing rows @@ -331,8 +382,42 @@ class FieldSpecificationImpl SortedArrayView< localIndex const > const & targetSet, arrayView1d< globalIndex const > const & dofMap, CRSMatrixView< real64, globalIndex const > const & matrix ); + +private: + + /** + * @brief Apply the lambda to each component of the field specification + * @tparam LAMBDA The type of lambda function passed into the parameter list. + * @param fs The field specification data object + * @param lambda The lambda being executed + */ + template< typename LAMBDA > + static void forEachComponent( FieldSpecification const & fs, LAMBDA && lambda ); + }; +template< typename LAMBDA > +void FieldSpecificationImpl::forEachComponent( FieldSpecification const & fs, LAMBDA && lambda ) +{ + if( fs.usesNonScalarValues() ) + { + localIndex const numComponents = fs.getScales().size(); + for( localIndex comp = 0; comp < numComponents; ++comp ) + { + // allow for a single functionName to be applied to every component + localIndex const fnIdx = ( fs.getFunctionName().size() == 1 ) ? localIndex( 0 ) : comp; + string const & functionName = fs.getFunctionNames().empty() ? string{} + : fs.getFunctionNames()[ fnIdx ]; + + lambda( comp, fs.getScales()[ comp ], functionName ); + } + } + else + { + lambda( fs.getComponent(), fs.getScale(), fs.getFunctionName() ); + } +} + template< typename OBJECT_TYPE, typename BC_TYPE, typename LAMBDA > void FieldSpecificationImpl::apply( BC_TYPE const & fs, MeshLevel & mesh, @@ -384,21 +469,24 @@ FieldSpecificationImpl:: real64 const time, dataRepository::Group & dataGroup ) { - integer const component = fs.getComponent(); - string const & functionName = fs.getFunctionName(); FunctionManager & functionManager = FunctionManager::getInstance(); - if( functionName.empty() ) + forEachComponent( fs, + [&]( integer const component, + real64 const scale, + string const & functionName ) { - real64 const value = fs.getScale(); - forAll< POLICY >( targetSet.size(), [=] GEOS_HOST_DEVICE ( localIndex const i ) + if( functionName.empty() ) { - localIndex const a = targetSet[ i ]; - FIELD_OP::SpecifyFieldValue( field, a, component, value ); - } ); - } - else - { + real64 const value = scale; + forAll< POLICY >( targetSet.size(), [=] GEOS_HOST_DEVICE ( localIndex const i ) + { + localIndex const a = targetSet[ i ]; + FIELD_OP::SpecifyFieldValue( field, a, component, value ); + } ); + return; + } + FunctionBase const & function = [&]() -> FunctionBase const & { try @@ -422,7 +510,7 @@ FieldSpecificationImpl:: if( function.isFunctionOfTime()==2 ) { - real64 const value = fs.getScale() * function.evaluate( &time ); + real64 const value = scale * function.evaluate( &time ); forAll< POLICY >( targetSet.size(), [=] GEOS_HOST_DEVICE ( localIndex const i ) { localIndex const a = targetSet[ i ]; @@ -434,14 +522,13 @@ FieldSpecificationImpl:: real64_array result( static_cast< localIndex >( targetSet.size() ) ); function.evaluate( dataGroup, time, targetSet, result ); arrayView1d< real64 const > const & resultView = result.toViewConst(); - real64 const scale = fs.getScale(); forAll< POLICY >( targetSet.size(), [=] GEOS_HOST_DEVICE ( localIndex const i ) { localIndex const a = targetSet[ i ]; FIELD_OP::SpecifyFieldValue( field, a, component, scale * resultView[i] ); } ); } - } + } ); } template< typename FIELD_OP, typename POLICY > @@ -486,20 +573,40 @@ FieldSpecificationImpl:: arrayView1d< real64 > const & rhs, ArrayView< T const, NDIM, USD > const & fieldView ) { - integer const component = fs.getComponent(); - applyBoundaryConditionToSystem< FIELD_OP, POLICY >( fs, - targetSet, - time, - dataGroup, - dofMap, - dofRankOffset, - matrix, - rhs, - [fieldView, component] GEOS_HOST_DEVICE ( localIndex const a ) + forEachComponent( fs, + [&]( integer const component, + real64 const scale, + string const & functionName ) { - real64 value = 0.0; - FieldSpecificationEqual::readFieldValue( fieldView, a, component, value ); - return value; + integer const comp = ( component >= 0 ) ? component : 0; + + array1d< globalIndex > dofArray( targetSet.size() ); + arrayView1d< globalIndex > const & dof = dofArray.toView(); + + array1d< real64 > rhsContributionArray( targetSet.size() ); + arrayView1d< real64 > const & rhsContribution = rhsContributionArray.toView(); + + computeRhsContribution< FIELD_OP, POLICY >( fs, + comp, + scale, + functionName, + targetSet, + time, + 1.0, + dataGroup, + dofMap, + dofRankOffset, + matrix, + dof, + rhsContribution, + [fieldView, component] GEOS_HOST_DEVICE ( localIndex const a ) + { + real64 value = 0.0; + FieldSpecificationEqual::readFieldValue( fieldView, a, component, value ); + return value; + } ); + + FIELD_OP::template prescribeRhsValues< POLICY >( rhs, dof, dofRankOffset, rhsContribution ); } ); } @@ -576,31 +683,45 @@ FieldSpecificationImpl:: arrayView1d< real64 > const & rhs, LAMBDA && lambda ) { - array1d< globalIndex > dofArray( targetSet.size() ); - arrayView1d< globalIndex > const & dof = dofArray.toView(); - - array1d< real64 > rhsContributionArray( targetSet.size() ); - arrayView1d< real64 > const & rhsContribution = rhsContributionArray.toView(); - - computeRhsContribution< FIELD_OP, POLICY, LAMBDA >( fs, - targetSet, - time, - dt, - dataGroup, - dofMap, - dofRankOffset, - matrix, - dof, - rhsContribution, - std::forward< LAMBDA >( lambda ) ); - - FIELD_OP::template prescribeRhsValues< POLICY >( rhs, dof, dofRankOffset, rhsContribution ); + forEachComponent( fs, + [&]( integer const component, + real64 const scale, + string const & functionName ) + { + integer const comp = ( component >= 0 ) ? component : 0; + + array1d< globalIndex > dofArray( targetSet.size() ); + arrayView1d< globalIndex > const & dof = dofArray.toView(); + + array1d< real64 > rhsContributionArray( targetSet.size() ); + arrayView1d< real64 > const & rhsContribution = rhsContributionArray.toView(); + + computeRhsContribution< FIELD_OP, POLICY >( fs, + comp, + scale, + functionName, + targetSet, + time, + dt, + dataGroup, + dofMap, + dofRankOffset, + matrix, + dof, + rhsContribution, + lambda ); + + FIELD_OP::template prescribeRhsValues< POLICY >( rhs, dof, dofRankOffset, rhsContribution ); + } ); } template< typename FIELD_OP, typename POLICY, typename LAMBDA > void FieldSpecificationImpl:: computeRhsContribution( FieldSpecification const & fs, + integer const component, + real64 const scale, + string const & functionName, SortedArrayView< localIndex const > const & targetSet, real64 const time, real64 const dt, @@ -612,8 +733,6 @@ FieldSpecificationImpl:: arrayView1d< real64 > const & rhsContribution, LAMBDA && lambda ) { - integer const component = ( fs.getComponent() >= 0 ) ? fs.getComponent() : 0; - string const & functionName = fs.getFunctionName(); FunctionManager & functionManager = FunctionManager::getInstance(); // Compute the value of the rhs terms, and collect the dof numbers @@ -623,7 +742,7 @@ FieldSpecificationImpl:: if( functionName.empty() || functionManager.getGroup< FunctionBase >( functionName ).isFunctionOfTime() == 2 ) { - real64 value = fs.getScale() * dt; + real64 value = scale * dt; if( !functionName.empty() ) { FunctionBase const & function = functionManager.getGroup< FunctionBase >( functionName ); @@ -658,7 +777,7 @@ FieldSpecificationImpl:: real64_array resultsArray( targetSet.size() ); function.evaluate( dataGroup, time, targetSet, resultsArray ); arrayView1d< real64 const > const & results = resultsArray.toViewConst(); - real64 const value = fs.getScale() * dt; + real64 const value = scale * dt; forAll< POLICY >( targetSet.size(), [targetSet, @@ -684,6 +803,39 @@ FieldSpecificationImpl:: } } +template< typename FIELD_OP, typename POLICY, typename LAMBDA > +void +FieldSpecificationImpl:: + computeRhsContribution( FieldSpecification const & fs, + SortedArrayView< localIndex const > const & targetSet, + real64 const time, + real64 const dt, + dataRepository::Group const & dataGroup, + arrayView1d< globalIndex const > const & dofMap, + globalIndex const dofRankOffset, + CRSMatrixView< real64, globalIndex const > const & matrix, + arrayView1d< globalIndex > const & dof, + arrayView1d< real64 > const & rhsContribution, + LAMBDA && lambda ) +{ + computeRhsContribution< FIELD_OP, POLICY, LAMBDA >( fs, + ( fs.getComponent() >= 0 ) ? fs.getComponent() : 0, + fs.getScale(), + fs.getFunctionName(), + targetSet, + time, + dt, + dataGroup, + dofMap, + dofRankOffset, + matrix, + dof, + rhsContribution, + std::forward< LAMBDA >( lambda ) ); + +} + + template< typename POLICY > void FieldSpecificationImpl:: @@ -692,20 +844,26 @@ FieldSpecificationImpl:: arrayView1d< globalIndex const > const & dofMap, CRSMatrixView< real64, globalIndex const > const & matrix ) { - integer const component = ( fs.getComponent() >= 0 ) ? fs.getComponent() : 0; - forAll< POLICY >( targetSet.size(), - [targetSet, dofMap, matrix, component] GEOS_HOST_DEVICE ( localIndex const i ) + forEachComponent( fs, + [&]( integer const component, + real64 const scale, + string const & functionName ) { - localIndex const a = targetSet[ i ]; - globalIndex const dof = dofMap[ a ] + component; + GEOS_UNUSED_VAR( scale, functionName ); + integer const comp = ( component >= 0 ) ? component : 0; + forAll< POLICY >( targetSet.size(), [targetSet, dofMap, matrix, comp] GEOS_HOST_DEVICE ( localIndex const i ) + { + localIndex const a = targetSet[ i ]; + globalIndex const dof = dofMap[ a ] + comp; - arraySlice1d< real64 > const entries = matrix.getEntries( dof ); - localIndex const numEntries = matrix.numNonZeros( dof ); + arraySlice1d< real64 > const entries = matrix.getEntries( dof ); + localIndex const numEntries = matrix.numNonZeros( dof ); - for( localIndex j = 0; j < numEntries; ++j ) - { - entries[ j ] = 0; - } + for( localIndex j = 0; j < numEntries; ++j ) + { + entries[ j ] = 0; + } + } ); } ); } diff --git a/src/coreComponents/fieldSpecification/PerfectlyMatchedLayer.cpp b/src/coreComponents/fieldSpecification/PerfectlyMatchedLayer.cpp index e8bfabcf5bd..765db5cb976 100644 --- a/src/coreComponents/fieldSpecification/PerfectlyMatchedLayer.cpp +++ b/src/coreComponents/fieldSpecification/PerfectlyMatchedLayer.cpp @@ -78,6 +78,8 @@ PerfectlyMatchedLayer::PerfectlyMatchedLayer( string const & name, Group * const void PerfectlyMatchedLayer::postInputInitialization() { + FieldSpecification::postInputInitialization(); + GEOS_THROW_IF( (m_xMax[0]( FieldSpecification::viewKeyStruct::functionNameString() ). - setDescription( GEOS_FMT( "Name of a function that specifies the variation of the production rate variations of this {}." + getWrapper< string_array >( FieldSpecification::viewKeyStruct::functionNameString() ). + setDescription( GEOS_FMT( "Name(s) of function(s) 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 produced fluid rate unit is in kg by default, or in mole if the flow solver uses moles.", catalogName(), FieldSpecification::viewKeyStruct::scaleString() ) ); - getWrapper< real64 >( FieldSpecification::viewKeyStruct::scaleString() ). - setDescription( GEOS_FMT( "Multiplier of the {0} value. If no {0} is provided, this value is used directly.", + getWrapper< array1d< real64 > >( FieldSpecification::viewKeyStruct::scaleString() ). + setDescription( GEOS_FMT( "Multiplier(s) of the {0} value(s). If no {0} is provided, value(s) are used directly.", FieldSpecification::viewKeyStruct::functionNameString() ) ); } diff --git a/src/coreComponents/fieldSpecification/TractionBoundaryCondition.cpp b/src/coreComponents/fieldSpecification/TractionBoundaryCondition.cpp index 6025988406e..30375f025c5 100644 --- a/src/coreComponents/fieldSpecification/TractionBoundaryCondition.cpp +++ b/src/coreComponents/fieldSpecification/TractionBoundaryCondition.cpp @@ -78,6 +78,8 @@ TractionBoundaryCondition::TractionBoundaryCondition( string const & name, Group void TractionBoundaryCondition::postInputInitialization() { + FieldSpecification::postInputInitialization(); + if( m_tractionType == TractionType::vector ) { GEOS_ERROR_IF( LvArray::tensorOps::l2Norm< 3 >( getDirection() ) < 1e-20, diff --git a/src/coreComponents/schema/schema.xsd b/src/coreComponents/schema/schema.xsd index 00589e978ab..b63b6d1e74c 100644 --- a/src/coreComponents/schema/schema.xsd +++ b/src/coreComponents/schema/schema.xsd @@ -30,7 +30,7 @@ - + @@ -55,7 +55,7 @@ - + @@ -65,7 +65,7 @@ - + @@ -85,7 +85,7 @@ - + @@ -110,7 +110,7 @@ - + @@ -120,7 +120,7 @@ - + @@ -140,7 +140,7 @@ - + @@ -165,7 +165,7 @@ - + @@ -1347,8 +1347,6 @@ This keyword is ignored for single-phase flow simulations--> - - @@ -1361,15 +1359,15 @@ When set to "warning", output a warning. When set to "error", output a throw. --> - - + + - - + + @@ -1383,8 +1381,6 @@ A set can be be defined by a 'Geometry' component, or correspond to imported set - - - - + + - - + + @@ -1419,8 +1415,6 @@ A set can be be defined by a 'Geometry' component, or correspond to imported set - - - - + + - - + + @@ -1455,8 +1449,6 @@ A set can be be defined by a 'Geometry' component, or correspond to imported set - - @@ -1481,8 +1473,8 @@ When set to "warning", output a warning. When set to "error", output a throw. --> - - + + @@ -1491,16 +1483,14 @@ When set to "error", output a throw. - - + + - - - - + + - - + + @@ -1544,8 +1534,6 @@ A set can be be defined by a 'Geometry' component, or correspond to imported set - - - - + + - - + + @@ -1577,8 +1565,6 @@ A set can be be defined by a 'Geometry' component, or correspond to imported set - - @@ -1591,8 +1577,8 @@ When set to "warning", output a warning. When set to "error", output a throw. --> - - + + @@ -1601,8 +1587,8 @@ When set to "error", output a throw. - - + +