feat: add way for field specification to scale input data

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

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();
 }
 

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;
+
 };
 
 

src/coreComponents/fieldSpecification/FieldSpecificationManager.cpp

@@ -216,7 +216,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,
@@ -225,10 +225,28 @@ void FieldSpecificationManager::applyInitialConditions( MeshLevel & mesh ) const
                                               Group & targetGroup,
                                               string const fieldName )
       {
-        bc.applyFieldValue< FieldSpecificationEqual >( targetSet, 0.0, targetGroup, fieldName );
+          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 */

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
@@ -251,6 +257,7 @@ FieldSpecificationManager::
               Group & targetGroup,
               string const & targetField )
   {
+
     fs.applyFieldValue< FieldSpecificationEqual, POLICY >( targetSet, time, targetGroup, targetField );
     lambda( fs, targetSet );
   } );

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 );
     } );