Merge pull request #13989 from KratosMultiphysics/optapp/vm_cload

sunethwarna · web-flow · commit 71fa255ccb97 · 2025-12-24T11:46:15.000-08:00
[OptApp] Adding cloud mesh option to explicit filtering
diff --git a/applications/OptimizationApplication/custom_python/add_custom_filters_to_python.cpp b/applications/OptimizationApplication/custom_python/add_custom_filters_to_python.cpp
@@ -39,7 +39,7 @@ void AddExplicitFilterUtils(
     namespace py = pybind11;
 
     py::class_<ExplicitFilterUtils<TContainerType>, typename ExplicitFilterUtils<TContainerType>::Pointer>(m, rName.c_str())
-        .def(py::init<const ModelPart&, const std::string&, const std::size_t, const std::size_t>(), py::arg("model_part"), py::arg("kernel_function_type"), py::arg("max_number_of_neighbours"), py::arg("echo_level"))
+        .def(py::init<const ModelPart&, const std::string&, const std::size_t, const std::size_t, bool>(), py::arg("model_part"), py::arg("kernel_function_type"), py::arg("max_number_of_neighbours"), py::arg("echo_level"), py::arg("node_cloud_mesh") = false)
         .def("SetRadius", &ExplicitFilterUtils<TContainerType>::SetRadius, py::arg("filter_radius"))
         .def("SetDamping", &ExplicitFilterUtils<TContainerType>::SetDamping, py::arg("damping"))
         .def("ForwardFilterField", &ExplicitFilterUtils<TContainerType>::ForwardFilterField, py::arg("mesh_independent_control_space_field"))
diff --git a/applications/OptimizationApplication/custom_utilities/filtering/explicit_filter_utils.cpp b/applications/OptimizationApplication/custom_utilities/filtering/explicit_filter_utils.cpp
@@ -65,25 +65,42 @@ Expression::ConstPointer GetNodalDomainSizeExpression(
     return p_expression;
 }
 
-template<class TEntityType>
-double GetDomainSize(
-    const EntityPoint<TEntityType>& rPoint,
-    Expression const * const pExpression)
+struct NodeCloudsType
 {
-    return rPoint.GetEntity().GetGeometry().DomainSize();
-}
+    template<class TEntityType>
+    static inline double GetDomainSize(
+        const EntityPoint<TEntityType>& rPoint,
+        Expression const * const pExpression)
+    {
+        if constexpr(!std::is_same_v<TEntityType, ModelPart::NodeType>) {
+            KRATOS_ERROR << "Point / node clouds should only be used with nodes.";
+        }
+        return 1.0;
+    }
 
-template<>
-double GetDomainSize(
-    const EntityPoint<ModelPart::NodeType>& rPoint,
-    Expression const * const pExpression)
-{
-    return pExpression->Evaluate(rPoint.Id(), rPoint.Id(), 0);
-}
+    static inline double Compute(
+        const double Value,
+        const double DomainSize)
+    {
+        return Value;
+    }
+};
 
 struct MeshIndependentType
 {
-    static double Compute(
+    template<class TEntityType>
+    static inline double GetDomainSize(
+        const EntityPoint<TEntityType>& rPoint,
+        Expression const * const pExpression)
+    {
+        if constexpr(std::is_same_v<TEntityType, ModelPart::NodeType>) {
+            return pExpression->Evaluate(rPoint.Id(), rPoint.Id(), 0);
+        } else {
+            return rPoint.GetEntity().GetGeometry().DomainSize();
+        }
+    }
+
+    static inline double Compute(
         const double Value,
         const double DomainSize)
     {
@@ -93,15 +110,27 @@ struct MeshIndependentType
 
 struct MeshDependentType
 {
-    static double Compute(
+    template<class TEntityType>
+    static inline double GetDomainSize(
+        const EntityPoint<TEntityType>& rPoint,
+        Expression const * const pExpression)
+    {
+        if constexpr(std::is_same_v<TEntityType, ModelPart::NodeType>) {
+            return pExpression->Evaluate(rPoint.Id(), rPoint.Id(), 0);
+        } else {
+            return rPoint.GetEntity().GetGeometry().DomainSize();
+        }
+    }
+
+    static inline double Compute(
         const double Value,
         const double DomainSize)
     {
         return Value / DomainSize;
     }
 };
 
-template <class TEntityType>
+template <class TMeshDependencyType, class TEntityType>
 void ComputeWeightForAllNeighbors(
     double& rSumOfWeights,
     std::vector<double>& rListOfWeights,
@@ -114,7 +143,7 @@ void ComputeWeightForAllNeighbors(
     Expression const * const pExpression)
 {
     for (IndexType i = 0; i < NumberOfNeighbours; ++i) {
-        const double domain_size = GetDomainSize(*rNeighbourNodes[i], pExpression);
+        const double domain_size = TMeshDependencyType::GetDomainSize(*rNeighbourNodes[i], pExpression);
         const double filter_weight = rFilterFunction.ComputeWeight(Radius, std::sqrt(rSquaredDistances[i])) * domain_size;
         rListOfWeights[i] = filter_weight;
         rSumOfWeights += filter_weight;
@@ -128,10 +157,12 @@ ExplicitFilterUtils<TContainerType>::ExplicitFilterUtils(
     const ModelPart& rModelPart,
     const std::string& rKernelFunctionType,
     const IndexType MaxNumberOfNeighbours,
-    const IndexType EchoLevel)
+    const IndexType EchoLevel,
+    const bool NodeCloudMesh)
     : mrModelPart(rModelPart),
       mMaxNumberOfNeighbors(MaxNumberOfNeighbours),
-      mEchoLevel(EchoLevel)
+      mEchoLevel(EchoLevel),
+      mNodeCloudMesh(NodeCloudMesh)
 {
     mpKernelFunction = Kratos::make_unique<FilterFunction>(rKernelFunctionType);
 }
@@ -184,20 +215,22 @@ void ExplicitFilterUtils<TContainerType>::ExplicitFilterUtils::Update()
     mpSearchTree =  Kratos::make_shared<ExplicitFilterUtils::KDTree>(mEntityPointVector.begin(), mEntityPointVector.end(), mBucketSize);
 
     if constexpr(std::is_same_v<TContainerType, ModelPart::NodesContainerType>) {
-        const auto& r_nodes = mrModelPart.Nodes();
-        const auto& r_conditions = mrModelPart.Conditions();
-        const IndexType number_of_conditions = r_conditions.size();
-
-        const auto& r_elements = mrModelPart.Elements();
-        const IndexType number_of_elements = r_elements.size();
-
-        if (number_of_elements > 0) {
-            mpNodalDomainSizeExpression = ExplicitFilterUtilsHelperUtilities::GetNodalDomainSizeExpression(r_elements, r_nodes);
-        } else if (number_of_conditions > 0) {
-            mpNodalDomainSizeExpression = ExplicitFilterUtilsHelperUtilities::GetNodalDomainSizeExpression(r_conditions, r_nodes);
-        } else {
-            KRATOS_ERROR << "Nodal mapping requires atleast either conditions or elements to be present in "
-                         << mrModelPart.FullName() << ".\n";
+        if (!mNodeCloudMesh) {
+            const auto& r_nodes = mrModelPart.Nodes();
+            const auto& r_conditions = mrModelPart.Conditions();
+            const IndexType number_of_conditions = r_conditions.size();
+
+            const auto& r_elements = mrModelPart.Elements();
+            const IndexType number_of_elements = r_elements.size();
+
+            if (number_of_elements > 0) {
+                mpNodalDomainSizeExpression = ExplicitFilterUtilsHelperUtilities::GetNodalDomainSizeExpression(r_elements, r_nodes);
+            } else if (number_of_conditions > 0) {
+                mpNodalDomainSizeExpression = ExplicitFilterUtilsHelperUtilities::GetNodalDomainSizeExpression(r_conditions, r_nodes);
+            } else {
+                KRATOS_ERROR << "Nodal mapping requires atleast either conditions or elements to be present in "
+                            << mrModelPart.FullName() << ".\n";
+            }
         }
     }
 
@@ -230,7 +263,8 @@ void ExplicitFilterUtils<TContainerType>::CheckField(const ContainerExpression<T
 }
 
 template<class TContainerType>
-ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::ForwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
+template<class TMeshDependencyType>
+ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::GenericForwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
 {
     KRATOS_TRY
 
@@ -268,7 +302,7 @@ ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::Forward
             << mMaxNumberOfNeighbors << " ].\n";
 
         double sum_of_weights = 0.0;
-        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors(
+        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors<TMeshDependencyType>(
             sum_of_weights, rTLS.mListOfWeights, *mpKernelFunction, radius,
             entity_point, rTLS.mNeighbourEntityPoints, rTLS.mResultingSquaredDistances, number_of_neighbors, this->mpNodalDomainSizeExpression.get());
 
@@ -293,6 +327,7 @@ ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::Forward
 
     KRATOS_CATCH("");
 }
+
 template<class TContainerType>
 template<class TMeshDependencyType>
 ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::GenericBackwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
@@ -337,14 +372,14 @@ ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::Generic
             << mMaxNumberOfNeighbors << " ].\n";
 
         double sum_of_weights = 0.0;
-        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors(
+        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors<TMeshDependencyType>(
             sum_of_weights, rTLS.mListOfWeights, *mpKernelFunction, radius,
             entity_point, rTLS.mNeighbourEntityPoints, rTLS.mResultingSquaredDistances, number_of_neighbors, this->mpNodalDomainSizeExpression.get());
 
         mpDamping->Apply(rTLS.mListOfDampedWeights, rTLS.mListOfWeights, Index, number_of_neighbors, rTLS.mNeighbourEntityPoints);
 
         const IndexType current_data_begin = Index * stride;
-        const double domain_size = ExplicitFilterUtilsHelperUtilities::GetDomainSize(entity_point, mpNodalDomainSizeExpression.get());
+        const double domain_size = TMeshDependencyType::GetDomainSize(entity_point, mpNodalDomainSizeExpression.get());
 
         for (IndexType j = 0; j < stride; ++j) {
             const auto& r_damped_weights = rTLS.mListOfDampedWeights[j];
@@ -369,19 +404,8 @@ ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::Generic
 }
 
 template<class TContainerType>
-ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::BackwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
-{
-    return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::MeshIndependentType>(rContainerExpression);
-}
-
-template<class TContainerType>
-ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::BackwardFilterIntegratedField(const ContainerExpression<TContainerType>& rContainerExpression) const
-{
-    return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::MeshDependentType>(rContainerExpression);
-}
-
-template<class TContainerType>
-void ExplicitFilterUtils<TContainerType>::GetIntegrationWeights(ContainerExpression<TContainerType>& rContainerExpression) const
+template<class TMeshDependencyType>
+void ExplicitFilterUtils<TContainerType>::GenericGetIntegrationWeights(ContainerExpression<TContainerType>& rContainerExpression) const
 {
 
     KRATOS_ERROR_IF_NOT(&rContainerExpression.GetModelPart() == &mrModelPart)
@@ -394,7 +418,7 @@ void ExplicitFilterUtils<TContainerType>::GetIntegrationWeights(ContainerExpress
 
     IndexPartition<IndexType>(r_container.size()).for_each([&](const IndexType Index){
         const EntityPoint<EntityType> entity(*(r_container.begin() + Index), Index);
-        const auto integration_weight = ExplicitFilterUtilsHelperUtilities::GetDomainSize(entity, this->mpNodalDomainSizeExpression.get());
+        const auto integration_weight = TMeshDependencyType::GetDomainSize(entity, this->mpNodalDomainSizeExpression.get());
         const IndexType current_data_begin = Index * stride;
         for (IndexType j = 0; j < stride; ++j) {
             double& current_index_value = *(p_expression->begin() + current_data_begin + j);
@@ -403,6 +427,46 @@ void ExplicitFilterUtils<TContainerType>::GetIntegrationWeights(ContainerExpress
     });
 }
 
+template<class TContainerType>
+ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::ForwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
+{
+    if (mNodeCloudMesh) {
+        return GenericForwardFilterField<ExplicitFilterUtilsHelperUtilities::NodeCloudsType>(rContainerExpression);
+    } else {
+        return GenericForwardFilterField<ExplicitFilterUtilsHelperUtilities::MeshIndependentType>(rContainerExpression);
+    }
+}
+
+template<class TContainerType>
+ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::BackwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const
+{
+    if (mNodeCloudMesh) {
+        return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::NodeCloudsType>(rContainerExpression);
+    } else {
+        return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::MeshIndependentType>(rContainerExpression);
+    }
+}
+
+template<class TContainerType>
+ContainerExpression<TContainerType> ExplicitFilterUtils<TContainerType>::BackwardFilterIntegratedField(const ContainerExpression<TContainerType>& rContainerExpression) const
+{
+    if (mNodeCloudMesh) {
+        return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::NodeCloudsType>(rContainerExpression);
+    } else {
+        return GenericBackwardFilterField<ExplicitFilterUtilsHelperUtilities::MeshDependentType>(rContainerExpression);
+    }
+}
+
+template<class TContainerType>
+void ExplicitFilterUtils<TContainerType>::GetIntegrationWeights(ContainerExpression<TContainerType>& rContainerExpression) const
+{
+    if (mNodeCloudMesh) {
+        GenericGetIntegrationWeights<ExplicitFilterUtilsHelperUtilities::NodeCloudsType>(rContainerExpression);
+    } else {
+        GenericGetIntegrationWeights<ExplicitFilterUtilsHelperUtilities::MeshIndependentType>(rContainerExpression);
+    }
+}
+
 template<class TContainerType>
 std::string ExplicitFilterUtils<TContainerType>::Info() const
 {
@@ -456,7 +520,7 @@ void ExplicitFilterUtils<TContainerType>::CalculateMatrix(Matrix& rOutput) const
 
         std::vector<double> list_of_weights(number_of_neighbors, 0.0);
         double sum_of_weights = 0.0;
-        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors(
+        ExplicitFilterUtilsHelperUtilities::ComputeWeightForAllNeighbors<ExplicitFilterUtilsHelperUtilities::MeshIndependentType>(
             sum_of_weights, list_of_weights, *mpKernelFunction, radius,
             *mEntityPointVector[Index], rTLS.mNeighbourEntityPoints, rTLS.mResultingSquaredDistances, number_of_neighbors, this->mpNodalDomainSizeExpression.get());
 
diff --git a/applications/OptimizationApplication/custom_utilities/filtering/explicit_filter_utils.h b/applications/OptimizationApplication/custom_utilities/filtering/explicit_filter_utils.h
@@ -63,7 +63,8 @@ class KRATOS_API(OPTIMIZATION_APPLICATION) ExplicitFilterUtils
         const ModelPart& rModelPart,
         const std::string& rKernelFunctionType,
         const IndexType MaxNumberOfNeighbours,
-        const IndexType EchoLevel);
+        const IndexType EchoLevel,
+        const bool NodeCloudMesh);
 
     ///@}
     ///@name Public operations
@@ -160,12 +161,20 @@ class KRATOS_API(OPTIMIZATION_APPLICATION) ExplicitFilterUtils
 
     typename KDTree::Pointer mpSearchTree;
 
+    bool mNodeCloudMesh;
+
     ///@}
     ///@name Private operations
     ///@{
 
     void CheckField(const ContainerExpression<TContainerType>& rContainerExpression) const;
 
+    template<class TMeshDependencyType>
+    void GenericGetIntegrationWeights(ContainerExpression<TContainerType>& rContainerExpression) const;
+
+    template<class TMeshDependencyType>
+    ContainerExpression<TContainerType> GenericForwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const;
+
     template<class TMeshDependencyType>
     ContainerExpression<TContainerType> GenericBackwardFilterField(const ContainerExpression<TContainerType>& rContainerExpression) const;
 
diff --git a/applications/OptimizationApplication/python_scripts/filtering/explicit_filter.py b/applications/OptimizationApplication/python_scripts/filtering/explicit_filter.py
@@ -17,6 +17,7 @@ class ExplicitFilter(Filter):
     def GetDefaultParameters(cls) -> Kratos.Parameters:
         return Kratos.Parameters("""{
             "filter_type"               : "explicit_filter",
+            "node_cloud_mesh"           : false,
             "filter_function_type"      : "linear",
             "max_nodes_in_filter_radius": 100000,
             "echo_level"                : 0,
@@ -67,11 +68,16 @@ def Initialize(self) -> None:
         filter_function_type = self.parameters["filter_function_type"].GetString()
         max_nodes_in_filter_radius = self.parameters["max_nodes_in_filter_radius"].GetInt()
         echo_level = self.parameters["echo_level"].GetInt()
+        node_cloud_mesh = self.parameters["node_cloud_mesh"].GetBool()
         if self.data_location in [Kratos.Globals.DataLocation.NodeHistorical, Kratos.Globals.DataLocation.NodeNonHistorical]:
-            self.filter_utils = KratosOA.NodeExplicitFilterUtils(self.model_part, filter_function_type, max_nodes_in_filter_radius, echo_level)
+            self.filter_utils = KratosOA.NodeExplicitFilterUtils(self.model_part, filter_function_type, max_nodes_in_filter_radius, echo_level, node_cloud_mesh)
         elif self.data_location == Kratos.Globals.DataLocation.Condition:
+            if node_cloud_mesh:
+                raise RuntimeError("\"node_cloud_mesh\" flag can be only used for filtering meshes with nodes.")
             self.filter_utils = KratosOA.ConditionExplicitFilterUtils(self.model_part, filter_function_type, max_nodes_in_filter_radius, echo_level)
         elif self.data_location == Kratos.Globals.DataLocation.Element:
+            if node_cloud_mesh:
+                raise RuntimeError("\"node_cloud_mesh\" flag can be only used for filtering meshes with nodes.")
             self.filter_utils = KratosOA.ElementExplicitFilterUtils(self.model_part, filter_function_type, max_nodes_in_filter_radius, echo_level)
 
         self.Update()
diff --git a/applications/OptimizationApplication/tests/filtering/explicit_filter_reference_sigmoidal_cloud_mesh.vtu b/applications/OptimizationApplication/tests/filtering/explicit_filter_reference_sigmoidal_cloud_mesh.vtu
diff --git a/applications/OptimizationApplication/tests/filtering/explicit_filters_tests.py b/applications/OptimizationApplication/tests/filtering/explicit_filters_tests.py
