Add the possibility to provide a characteristic function

The characteristic function defines the domain of vertices that can be displaced.

src/language/modules/SchemeModule.cpp

@@ -131,6 +131,46 @@ SchemeModule::SchemeModule()
 
                               ));
 
+  this->_addBuiltinFunction("randomizeMesh",
+                            std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
+                              const IMesh>(std::shared_ptr<const IMesh>,
+                                           const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
+                                           const FunctionSymbolId&)>>(
+
+                              [](std::shared_ptr<const IMesh> p_mesh,
+                                 const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
+                                   bc_descriptor_list,
+                                 const FunctionSymbolId& function_symbol_id) -> std::shared_ptr<const IMesh> {
+                                switch (p_mesh->dimension()) {
+                                case 1: {
+                                  constexpr size_t Dimension = 1;
+                                  using MeshType             = Mesh<Connectivity<Dimension>>;
+                                  const MeshType& mesh       = dynamic_cast<const MeshType&>(*p_mesh);
+                                  MeshRandomizer randomizer(mesh, bc_descriptor_list);
+                                  return randomizer.getRandomizedMesh(function_symbol_id);
+                                }
+                                case 2: {
+                                  constexpr size_t Dimension = 2;
+                                  using MeshType             = Mesh<Connectivity<Dimension>>;
+                                  const MeshType& mesh       = dynamic_cast<const MeshType&>(*p_mesh);
+                                  MeshRandomizer randomizer(mesh, bc_descriptor_list);
+                                  return randomizer.getRandomizedMesh(function_symbol_id);
+                                }
+                                case 3: {
+                                  constexpr size_t Dimension = 3;
+                                  using MeshType             = Mesh<Connectivity<Dimension>>;
+                                  const MeshType& mesh       = dynamic_cast<const MeshType&>(*p_mesh);
+                                  MeshRandomizer randomizer(mesh, bc_descriptor_list);
+                                  return randomizer.getRandomizedMesh(function_symbol_id);
+                                }
+                                default: {
+                                  throw UnexpectedError("invalid mesh dimension");
+                                }
+                                }
+                              }
+
+                              ));
+
   this->_addBuiltinFunction("boundaryName",
                             std::make_shared<
                               BuiltinFunctionEmbedder<std::shared_ptr<const IBoundaryDescriptor>(const std::string&)>>(

src/mesh/MeshRandomizer.hpp

@@ -3,6 +3,7 @@
 
 #include <algebra/TinyMatrix.hpp>
 #include <algebra/TinyVector.hpp>
+#include <language/utils/InterpolateItemValue.hpp>
 #include <mesh/Connectivity.hpp>
 #include <mesh/ItemValueUtils.hpp>
 #include <mesh/Mesh.hpp>
@@ -137,9 +138,8 @@ class MeshRandomizer
     }
   }
 
- public:
-  std::shared_ptr<const MeshType>
-  getRandomizedMesh() const
+  NodeValue<Rd>
+  _getDisplacement() const
   {
     const ConnectivityType& connectivity = m_given_mesh.connectivity();
     NodeValue<const Rd> given_xr         = m_given_mesh.xr();
@@ -172,7 +172,7 @@ class MeshRandomizer
 
     synchronize(max_delta_xr);
 
-    std::uniform_real_distribution<> dis(-0.3, 0.3);
+    std::uniform_real_distribution<> dis(-0.5, 0.5);
 
     NodeValue<const int> node_numbers = connectivity.nodeNumber();
     using IdCorrespondance            = std::pair<int, NodeId>;
@@ -187,7 +187,7 @@ class MeshRandomizer
 
     RandomEngine& random_engine = RandomEngine::instance();
 
-    NodeValue<Rd> xr{connectivity};
+    NodeValue<Rd> shift_r{connectivity};
 
     int i_node_number = 0;
     for (size_t i = 0; i < node_numbers_to_node_id.size(); ++i) {
@@ -205,7 +205,7 @@ class MeshRandomizer
         shift[i_component] = max_delta * dis(random_engine.engine());
       }
 
-      xr[node_id] = shift;
+      shift_r[node_id] = shift;
 
       ++i_node_number;
     }
@@ -220,10 +220,37 @@ class MeshRandomizer
       }
     }
 
-    this->_applyBC(xr);
+    this->_applyBC(shift_r);
+    return shift_r;
+  }
+
+ public:
+  std::shared_ptr<const MeshType>
+  getRandomizedMesh() const
+  {
+    NodeValue<const Rd> given_xr = m_given_mesh.xr();
+
+    NodeValue<Rd> xr = this->_getDisplacement();
+
+    parallel_for(
+      m_given_mesh.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { xr[node_id] += given_xr[node_id]; });
+
+    return std::make_shared<MeshType>(m_given_mesh.shared_connectivity(), xr);
+  }
+
+  std::shared_ptr<const MeshType>
+  getRandomizedMesh(const FunctionSymbolId& function_symbol_id) const
+  {
+    NodeValue<const Rd> given_xr = m_given_mesh.xr();
+
+    NodeValue<const bool> is_displaced =
+      InterpolateItemValue<bool(const Rd)>::interpolate(function_symbol_id, given_xr);
+
+    NodeValue<Rd> xr = this->_getDisplacement();
 
     parallel_for(
-      connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) { xr[node_id] += given_xr[node_id]; });
+      m_given_mesh.numberOfNodes(),
+      PUGS_LAMBDA(const NodeId node_id) { xr[node_id] = is_displaced[node_id] * xr[node_id] + given_xr[node_id]; });
 
     return std::make_shared<MeshType>(m_given_mesh.shared_connectivity(), xr);
   }