From 18b8899a4202e6d80e6e328864b4cbbc134f539f Mon Sep 17 00:00:00 2001 From: Ryan Aronson Date: Mon, 2 Dec 2024 16:52:35 +0000 Subject: [PATCH 1/5] initial try --- .../common/FieldSpecificationOps.hpp | 84 +++++++++++++++++++ .../FieldSpecificationBase.cpp | 5 ++ .../FieldSpecificationBase.hpp | 14 ++++ .../FieldSpecificationManager.hpp | 16 +++- 4 files changed, 117 insertions(+), 2 deletions(-) diff --git a/src/coreComponents/common/FieldSpecificationOps.hpp b/src/coreComponents/common/FieldSpecificationOps.hpp index 534768330fd..287c165736b 100644 --- a/src/coreComponents/common/FieldSpecificationOps.hpp +++ b/src/coreComponents/common/FieldSpecificationOps.hpp @@ -67,6 +67,26 @@ struct OpAdd } }; +/** + * @brief OpMultiply Operator that multiplies by a value + */ +struct OpMultiply +{ + /** + * @brief Pointwise update of a value + * @tparam T type of the left-hand side + * @tparam U type of the right-hand side + * @param[in] lhs value to update + * @param[in] rhs input value + */ + template< typename T, typename U > + GEOS_HOST_DEVICE static inline + void apply( T & lhs, U const & rhs ) + { + lhs *= static_cast< T >( rhs ); + } +}; + /** * @brief FieldSpecificationOp */ @@ -681,6 +701,70 @@ struct FieldSpecificationAdd : public FieldSpecificationOp< OpAdd > }; + +/** + * @struct FieldSpecificationMultiply + * this struct a collection of static functions which adhere to an assumed interface for multiplying + * a value for a field. + */ +struct FieldSpecificationMultiply : public FieldSpecificationOp< OpMultiply > +{ + /// Alias for FieldSpecificationOp< OpMultiply > + using base_type = FieldSpecificationOp< OpMultiply >; + using base_type::SpecifyFieldValue; + + /** + * @brief Function to apply a value to a vector field for a single dof. + * @param[in] dof The degree of freedom that is to be modified. + * @param[in] dofRankOffset offset of dof indices on current rank + * @param[in] matrix A ParalleMatrix object: the system matrix. + * @param[out] rhs The rhs contribution to be modified + * @param[in] bcValue The value to multiply to rhs + * @param[in] fieldValue unused. + * + */ + GEOS_HOST_DEVICE + static inline void + SpecifyFieldValue( globalIndex const dof, + globalIndex const dofRankOffset, + CRSMatrixView< real64, globalIndex const > const & matrix, + real64 & rhs, + real64 const bcValue, + real64 const fieldValue ) + { + GEOS_UNUSED_VAR( dof ); + GEOS_UNUSED_VAR( dofRankOffset ); + GEOS_UNUSED_VAR( matrix ); + GEOS_UNUSED_VAR( fieldValue ); + rhs *= bcValue; + } + + /** + * @brief Function to add some values of a vector. + * @tparam POLICY the execution policy to use when setting values + * @param rhs the target right-hand side vector + * @param dof a list of global DOF indices to be set + * @param dofRankOffset offset of dof indices on current rank + * @param values a list of values corresponding to \p dof that will be added to \p rhs. + */ + template< typename POLICY > + static inline void prescribeRhsValues( arrayView1d< real64 > const & rhs, + arrayView1d< globalIndex const > const & dof, + globalIndex const dofRankOffset, + arrayView1d< real64 const > const & values ) + { + GEOS_ASSERT_EQ( dof.size(), values.size() ); + forAll< POLICY >( dof.size(), [rhs, dof, dofRankOffset, values] GEOS_HOST_DEVICE ( localIndex const a ) + { + globalIndex const localRow = dof[ a ] - dofRankOffset; + if( localRow >= 0 && localRow < rhs.size() ) + { rhs[ localRow ] *= values[ a ]; } + } ); + } + +}; + + } //namespace geos #endif //GEOS_COMMON_FIELDSPECIFICATIONOPS_HPP diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationBase.cpp b/src/coreComponents/fieldSpecification/FieldSpecificationBase.cpp index 45e58177dd4..7eaa924cc4e 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationBase.cpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationBase.cpp @@ -83,6 +83,11 @@ FieldSpecificationBase::FieldSpecificationBase( string const & name, Group * par setInputFlag( InputFlags::OPTIONAL ). setDescription( "Time at which the boundary condition will stop being applied." ); + registerWrapper( viewKeyStruct::isScalingString(), &m_isScaling ). + setApplyDefaultValue( 0 ). + setInputFlag( InputFlags::OPTIONAL ). + setDescription( "Boundary condition is a multiplicative scaling of the values already present." ); + enableLogLevelInput(); } diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationBase.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationBase.hpp index 8306eafca0b..ce67d21b33a 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationBase.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationBase.hpp @@ -392,6 +392,8 @@ class FieldSpecificationBase : public dataRepository::Group constexpr static char const * beginTimeString() { return "beginTime"; } /// @return The key for endTime constexpr static char const * endTimeString() { return "endTime"; } + /// @return The key for isScaling + constexpr static char const * isScalingString() { return "isScaling"; } }; /** @@ -476,6 +478,15 @@ class FieldSpecificationBase : public dataRepository::Group return m_initialCondition; } + /** + * Accessor + * @return const m_isScaling + */ + int isScaling() const + { + return m_isScaling; + } + /** * Accessor * @return const m_scale @@ -591,6 +602,9 @@ class FieldSpecificationBase : public dataRepository::Group /// The name of a function used to turn on and off the boundary condition. string m_bcApplicationFunctionName; + /// Whether or not the boundary condition is a multiplicative scaling of what is already present (from mesh) + int m_isScaling; + }; diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp index 91672902314..c8de57587af 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp @@ -251,7 +251,13 @@ FieldSpecificationManager:: Group & targetGroup, string const & targetField ) { - fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); + if (fs.isScaling()) + { + fs.applyFieldValue< FieldSpecificationMultiply, POLICY >( targetSet, time, targetGroup, targetField ); + } + else{ + fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); + } lambda( fs, targetSet ); } ); } @@ -275,7 +281,13 @@ FieldSpecificationManager:: string const & targetField ) { preLambda( fs, targetSet ); - fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); + if (fs.isScaling()) + { + fs.applyFieldValue< FieldSpecificationMultiply, POLICY >( targetSet, time, targetGroup, targetField ); + } + else{ + fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); + } postLambda( fs, targetSet ); } ); } From 6a260f7bd4f537c71e39eaa343687a31d955678e Mon Sep 17 00:00:00 2001 From: Ryan Aronson Date: Mon, 2 Dec 2024 11:19:57 -0800 Subject: [PATCH 2/5] missed the initial condition application, but now it is doing the mult 3 times --- .../fieldSpecification/FieldSpecificationManager.cpp | 9 ++++++++- .../fieldSpecification/FieldSpecificationManager.hpp | 4 +++- 2 files changed, 11 insertions(+), 2 deletions(-) diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp index 6edc963b457..00fc7950b41 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp @@ -224,7 +224,14 @@ void FieldSpecificationManager::applyInitialConditions( MeshLevel & mesh ) const Group & targetGroup, string const fieldName ) { - bc.applyFieldValue< FieldSpecificationEqual >( targetSet, 0.0, targetGroup, fieldName ); + if (fs.isScaling()) + { + bc.applyFieldValue< FieldSpecificationMultiply >( targetSet, 0.0, targetGroup, fieldName ); + } + else + { + bc.applyFieldValue< FieldSpecificationEqual >( targetSet, 0.0, targetGroup, fieldName ); + } } ); } } ); diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp index c8de57587af..09c5fddf84d 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp @@ -251,11 +251,13 @@ FieldSpecificationManager:: Group & targetGroup, string const & targetField ) { + if (fs.isScaling()) { fs.applyFieldValue< FieldSpecificationMultiply, POLICY >( targetSet, time, targetGroup, targetField ); } - else{ + else + { fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); } lambda( fs, targetSet ); From 61db5ac54d3f2df1c2a4ee8c6aa5f6338dbfe66a Mon Sep 17 00:00:00 2001 From: Ryan Aronson Date: Thu, 5 Dec 2024 16:53:57 -0800 Subject: [PATCH 3/5] fix multiple multiplications by only applying to one mesh level. works for internal mesh --- .../FieldSpecificationManager.cpp | 27 +++++++++++++------ .../FieldSpecificationManager.hpp | 6 +++++ .../mainInterface/ProblemManager.cpp | 1 + 3 files changed, 26 insertions(+), 8 deletions(-) diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp index 00fc7950b41..7103ed7cfe9 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp @@ -215,7 +215,7 @@ void FieldSpecificationManager::applyInitialConditions( MeshLevel & mesh ) const { this->forSubGroups< FieldSpecificationBase >( [&] ( FieldSpecificationBase const & fs ) { - if( fs.initialCondition() ) + if( fs.initialCondition() && !fs.isScaling() ) { fs.apply< dataRepository::Group >( mesh, [&]( FieldSpecificationBase const & bc, @@ -224,17 +224,28 @@ void FieldSpecificationManager::applyInitialConditions( MeshLevel & mesh ) const Group & targetGroup, string const fieldName ) { - if (fs.isScaling()) - { - bc.applyFieldValue< FieldSpecificationMultiply >( targetSet, 0.0, targetGroup, fieldName ); - } - else - { bc.applyFieldValue< FieldSpecificationEqual >( targetSet, 0.0, targetGroup, fieldName ); - } } ); } } ); } +void FieldSpecificationManager::applyScalingInitialConditions( MeshLevel & mesh ) const +{ + this->forSubGroups< FieldSpecificationBase >( [&] ( FieldSpecificationBase const & fs ) + { + if( fs.initialCondition() && fs.isScaling() ) + { + fs.apply< dataRepository::Group >( mesh, + [&]( FieldSpecificationBase const & bc, + string const &, + SortedArrayView< localIndex const > const & targetSet, + Group & targetGroup, + string const fieldName ) + { + bc.applyFieldValue< FieldSpecificationMultiply >( targetSet, 0.0, targetGroup, fieldName ); + } ); + } + } ); +} } /* namespace geos */ diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp index 09c5fddf84d..769a50660c4 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp @@ -178,6 +178,12 @@ class FieldSpecificationManager : public dataRepository::Group */ void applyInitialConditions( MeshLevel & mesh ) const; + /** + * @brief function to apply initial conditions which involve scaling a mesh field + * @param mesh the MeshLevel object + */ + void applyScalingInitialConditions( MeshLevel & mesh ) const; + /** * @brief function to validate the application of boundary conditions * @param mesh the MeshLevel object diff --git a/src/coreComponents/mainInterface/ProblemManager.cpp b/src/coreComponents/mainInterface/ProblemManager.cpp index d60a61663b7..d3e77dc2d18 100644 --- a/src/coreComponents/mainInterface/ProblemManager.cpp +++ b/src/coreComponents/mainInterface/ProblemManager.cpp @@ -1159,6 +1159,7 @@ void ProblemManager::applyInitialConditions() if( !meshLevel.isShallowCopy() ) // to avoid messages printed three times { m_fieldSpecificationManager->validateBoundaryConditions( meshLevel ); + m_fieldSpecificationManager->applyScalingInitialConditions( meshLevel );// with current implementation, scaling only works with shallow copy mesh levels } m_fieldSpecificationManager->applyInitialConditions( meshLevel ); } ); From a51495f2e1a3df201ca238703382694d05a2140f Mon Sep 17 00:00:00 2001 From: Ryan Aronson Date: Thu, 5 Dec 2024 17:09:59 -0800 Subject: [PATCH 4/5] Update FieldSpecificationManager.hpp --- .../FieldSpecificationManager.hpp | 19 +++---------------- 1 file changed, 3 insertions(+), 16 deletions(-) diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp index 769a50660c4..cd96cc743b0 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp @@ -258,14 +258,7 @@ FieldSpecificationManager:: string const & targetField ) { - if (fs.isScaling()) - { - fs.applyFieldValue< FieldSpecificationMultiply, POLICY >( targetSet, time, targetGroup, targetField ); - } - else - { - fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); - } + fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); lambda( fs, targetSet ); } ); } @@ -289,14 +282,8 @@ FieldSpecificationManager:: string const & targetField ) { preLambda( fs, targetSet ); - if (fs.isScaling()) - { - fs.applyFieldValue< FieldSpecificationMultiply, POLICY >( targetSet, time, targetGroup, targetField ); - } - else{ - fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); - } - postLambda( fs, targetSet ); + fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); + } ); } From 6c421b1ad2ee94426843bdd5de3c3a1ad5359aea Mon Sep 17 00:00:00 2001 From: Ryan Aronson Date: Thu, 5 Dec 2024 17:10:45 -0800 Subject: [PATCH 5/5] Update FieldSpecificationManager.hpp --- .../fieldSpecification/FieldSpecificationManager.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp index cd96cc743b0..c4832580629 100644 --- a/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp +++ b/src/coreComponents/fieldSpecification/FieldSpecificationManager.hpp @@ -283,7 +283,7 @@ FieldSpecificationManager:: { preLambda( fs, targetSet ); fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField ); - + postLambda( fs, targetSet ); } ); }