From 8e7381014c0bddc81b5eb2676e26a52638cd3079 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <github@naumov.de>
Date: Mon, 22 Mar 2021 11:38:50 +0100
Subject: [PATCH] [PL/DS] Use active elements to construct BC.

The set of deactivated nodes is now constructed on the fly using the
active elements set passed from the process variable.
---
 .../DeactivatedSubdomainDirichlet.cpp         | 26 +++++++++++++++++--
 .../DeactivatedSubdomainDirichlet.h           |  2 ++
 ProcessLib/ProcessVariable.cpp                |  1 +
 3 files changed, 27 insertions(+), 2 deletions(-)

diff --git a/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.cpp b/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.cpp
index 1fcf4942d99..3a6a0e20b23 100644
--- a/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.cpp
+++ b/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.cpp
@@ -20,6 +20,7 @@
 namespace ProcessLib
 {
 DeactivatedSubdomainDirichlet::DeactivatedSubdomainDirichlet(
+    std::vector<std::size_t> const* const active_element_ids,
     BaseLib::TimeInterval const& time_interval,
     ParameterLib::Parameter<double> const& parameter,
     DeactivatedSubdomainMesh const& subdomain,
@@ -29,7 +30,8 @@ DeactivatedSubdomainDirichlet::DeactivatedSubdomainDirichlet(
       _subdomain(subdomain),
       _variable_id(variable_id),
       _component_id(component_id),
-      _time_interval(time_interval)
+      _time_interval(time_interval),
+      _active_element_ids(active_element_ids)
 {
     config(dof_table_bulk);
 }
@@ -53,10 +55,30 @@ void DeactivatedSubdomainDirichlet::getEssentialBCValues(
     const double t, GlobalVector const& x,
     NumLib::IndexValueVector<GlobalIndexType>& bc_values) const
 {
+    auto const& bulk_element_ids =
+        *_subdomain.mesh->getProperties()
+             .template getPropertyVector<std::size_t>(
+                 "bulk_element_ids", MeshLib::MeshItemType::Cell, 1);
+
+    auto is_inactive = [&](MeshLib::Element const* const e) {
+        return none_of(
+            begin(*_active_element_ids), end(*_active_element_ids),
+            [&](auto const id) { return id == bulk_element_ids[e->getID()]; });
+    };
+
+    std::vector<MeshLib::Node*> inactive_nodes_in_bc_mesh;
+    std::copy_if(begin(_subdomain.inner_nodes), end(_subdomain.inner_nodes),
+                 back_inserter(inactive_nodes_in_bc_mesh),
+                 [&](MeshLib::Node* const n) {
+                     const auto& connected_elements = n->getElements();
+
+                     return std::all_of(begin(connected_elements),
+                                        end(connected_elements), is_inactive);
+                 });
     if (_time_interval.contains(t))
     {
         getEssentialBCValuesLocal(_parameter, *_subdomain.mesh,
-                                  _subdomain.inner_nodes, *_dof_table_boundary,
+                                  inactive_nodes_in_bc_mesh, *_dof_table_boundary,
                                   _variable_id, _component_id, t, x, bc_values);
         return;
     }
diff --git a/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.h b/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.h
index 7166e05e56f..cadfba543ad 100644
--- a/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.h
+++ b/ProcessLib/BoundaryCondition/DeactivatedSubdomainDirichlet.h
@@ -41,6 +41,7 @@ class DeactivatedSubdomainDirichlet final : public BoundaryCondition
 {
 public:
     DeactivatedSubdomainDirichlet(
+        std::vector<std::size_t> const* active_element_ids,
         BaseLib::TimeInterval const& time_interval,
         ParameterLib::Parameter<double> const& parameter,
         DeactivatedSubdomainMesh const& subdomain,
@@ -64,5 +65,6 @@ private:
     int const _component_id;
 
     BaseLib::TimeInterval const _time_interval;
+    std::vector<std::size_t> const* _active_element_ids = nullptr;
 };
 }  // namespace ProcessLib
diff --git a/ProcessLib/ProcessVariable.cpp b/ProcessLib/ProcessVariable.cpp
index 9d7e8b53a21..310d9937d4f 100644
--- a/ProcessLib/ProcessVariable.cpp
+++ b/ProcessLib/ProcessVariable.cpp
@@ -253,6 +253,7 @@ void ProcessVariable::createBoundaryConditionsForDeactivatedSubDomains(
                  component_id++)
             {
                 auto bc = std::make_unique<DeactivatedSubdomainDirichlet>(
+                    &_ids_of_active_elements,
                     deactivated_subdomain->time_interval, parameter,
                     *deactivated_subdomain_mesh, dof_table, variable_id,
                     component_id);
-- 
GitLab