diff --git a/ProcessLib/BoundaryCondition/BoundaryCondition.cpp b/ProcessLib/BoundaryCondition/BoundaryCondition.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..29d94a7fa0eb214af2aae30554092b532611ca81
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/BoundaryCondition.cpp
@@ -0,0 +1,46 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "BoundaryCondition.h"
+#include "MeshGeoToolsLib/BoundaryElementsSearcher.h"
+#include "MeshGeoToolsLib/MeshNodeSearcher.h"
+#include "BoundaryConditionConfig.h"
+#include "UniformDirichletBoundaryCondition.h"
+#include "UniformNeumannBoundaryCondition.h"
+
+namespace ProcessLib
+{
+std::unique_ptr<BoundaryCondition> createBoundaryCondition(
+ const BoundaryConditionConfig& config,
+ const NumLib::LocalToGlobalIndexMap& dof_table,
+ const MeshLib::Mesh& mesh,
+ const int variable_id,
+ const unsigned integration_order)
+{
+ MeshGeoToolsLib::MeshNodeSearcher& mesh_node_searcher =
+ MeshGeoToolsLib::MeshNodeSearcher::getMeshNodeSearcher(mesh);
+
+ MeshGeoToolsLib::BoundaryElementsSearcher boundary_element_searcher(
+ mesh, mesh_node_searcher);
+
+ auto const type = config.config.peekConfigParameter<std::string>("type");
+
+ if (type == "UniformDirichlet") {
+ return createUniformDirichletBoundaryCondition(
+ config, mesh_node_searcher, dof_table, variable_id);
+ } else if (type == "UniformNeumann") {
+ return createUniformNeumannBoundaryCondition(
+ config, boundary_element_searcher, dof_table, variable_id,
+ integration_order, mesh.getDimension());
+ } else {
+ OGS_FATAL("Unknown boundary condition type: `%s'.", type.c_str());
+ }
+}
+
+} // ProcessLib
diff --git a/ProcessLib/BoundaryCondition/BoundaryCondition.h b/ProcessLib/BoundaryCondition/BoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..674bcd324840377b3ac6b43f8e73f551cca73cc9
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/BoundaryCondition.h
@@ -0,0 +1,46 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_BOUNDARYCONDITION_H
+#define PROCESSLIB_BOUNDARYCONDITION_H
+
+#include "NumLib/NumericsConfig.h"
+
+namespace MeshLib
+{
+class Mesh;
+}
+
+namespace NumLib
+{
+class LocalToGlobalIndexMap;
+}
+
+namespace ProcessLib
+{
+struct BoundaryConditionConfig;
+
+class BoundaryCondition
+{
+public:
+ virtual void apply(const double t, GlobalVector const& x, GlobalMatrix& K,
+ GlobalVector& b) = 0;
+ virtual ~BoundaryCondition() = default;
+};
+
+std::unique_ptr<BoundaryCondition> createBoundaryCondition(
+ const BoundaryConditionConfig& config,
+ const NumLib::LocalToGlobalIndexMap& dof_table,
+ const MeshLib::Mesh& mesh,
+ const int variable_id,
+ const unsigned integration_order);
+
+} // ProcessLib
+
+#endif // PROCESSLIB_BOUNDARYCONDITION_H
diff --git a/ProcessLib/BoundaryCondition/BoundaryConditionCollection.cpp b/ProcessLib/BoundaryCondition/BoundaryConditionCollection.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..759e09adc0d222af886665e5498ecd62573f6241
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/BoundaryConditionCollection.cpp
@@ -0,0 +1,57 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "BoundaryConditionCollection.h"
+
+void initializeDirichletBCs(
+ std::vector<std::unique_ptr<ProcessLib::BoundaryCondition>> const&
+ boundary_conditions,
+ std::vector<NumLib::IndexValueVector<GlobalIndexType>>&
+ dirichlet_bcs)
+{
+ for (auto const& bc : boundary_conditions) {
+ if (auto* dirichlet_bc =
+ dynamic_cast<ProcessLib::DirichletBoundaryCondition*>(
+ bc.get())) {
+ dirichlet_bcs.emplace_back(dirichlet_bc->getBCValues());
+ }
+ }
+}
+
+namespace ProcessLib
+{
+void BoundaryConditionCollection::apply(const double t, GlobalVector const& x,
+ GlobalMatrix& K,
+ GlobalVector& b)
+{
+ for (auto const& bc : _boundary_conditions)
+ bc->apply(t, x, K, b);
+}
+
+void BoundaryConditionCollection::addBCsForProcessVariables(
+ std::vector<std::reference_wrapper<ProcessVariable>> const&
+ process_variables,
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ unsigned const integration_order)
+{
+ for (int variable_id = 0;
+ variable_id < static_cast<int>(process_variables.size());
+ ++variable_id)
+ {
+ ProcessVariable& pv = process_variables[variable_id];
+ auto bcs = pv.getBoundaryConditions(dof_table, variable_id,
+ integration_order);
+
+ std::move(bcs.begin(), bcs.end(),
+ std::back_inserter(_boundary_conditions));
+ }
+
+ initializeDirichletBCs(_boundary_conditions, _dirichlet_bcs);
+}
+}
diff --git a/ProcessLib/BoundaryCondition/BoundaryConditionCollection.h b/ProcessLib/BoundaryCondition/BoundaryConditionCollection.h
new file mode 100644
index 0000000000000000000000000000000000000000..48e899c594cc829af7de26f4b51dcf3a99c24392
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/BoundaryConditionCollection.h
@@ -0,0 +1,46 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_BOUNDARYCONDITIONCOLLECTION_H
+#define PROCESSLIB_BOUNDARYCONDITIONCOLLECTION_H
+
+#include "DirichletBoundaryCondition.h"
+#include "ProcessLib/ProcessVariable.h"
+
+namespace ProcessLib
+{
+
+class BoundaryConditionCollection : public BoundaryCondition
+{
+public:
+ void apply(const double t, GlobalVector const& x, GlobalMatrix& K,
+ GlobalVector& b) override;
+
+ std::vector<NumLib::IndexValueVector<GlobalIndexType>> const*
+ getKnownSolutions(double const /*t*/) const
+ {
+ // TODO time-dependent Dirichlet BCs.
+ return &_dirichlet_bcs;
+ }
+
+ void addBCsForProcessVariables(
+ std::vector<std::reference_wrapper<ProcessVariable>> const&
+ process_variables,
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ unsigned const integration_order);
+
+private:
+ std::vector<NumLib::IndexValueVector<GlobalIndexType>> _dirichlet_bcs;
+ std::vector<std::unique_ptr<BoundaryCondition>> _boundary_conditions;
+};
+
+
+} // ProcessLib
+
+#endif // PROCESSLIB_BOUNDARYCONDITIONCOLLECTION_H
diff --git a/ProcessLib/BoundaryCondition/BoundaryConditionConfig.h b/ProcessLib/BoundaryCondition/BoundaryConditionConfig.h
new file mode 100644
index 0000000000000000000000000000000000000000..b9c6d98e350dc747939bcca89fb67576d461e0de
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/BoundaryConditionConfig.h
@@ -0,0 +1,44 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_BOUNDARYCONDITIONCONFIG_H
+#define PROCESSLIB_BOUNDARYCONDITIONCONFIG_H
+
+#include "BaseLib/ConfigTree.h"
+#include "GeoLib/GEOObjects.h"
+
+namespace ProcessLib
+{
+
+struct BoundaryConditionConfig final
+{
+ BoundaryConditionConfig(BaseLib::ConfigTree&& config_,
+ GeoLib::GeoObject const& geometry_,
+ int const component_id_)
+ : config(std::move(config_)),
+ geometry(geometry_),
+ component_id(component_id_)
+ {
+ }
+
+ BoundaryConditionConfig(BoundaryConditionConfig&& other)
+ : config(std::move(other.config)),
+ geometry(other.geometry),
+ component_id(other.component_id)
+ {
+ }
+
+ BaseLib::ConfigTree config;
+ GeoLib::GeoObject const& geometry;
+ int const component_id;
+};
+
+} // ProcessLib
+
+#endif // PROCESSLIB_BOUNDARYCONDITIONCONFIG_H
diff --git a/ProcessLib/BoundaryCondition/CMakeLists.txt b/ProcessLib/BoundaryCondition/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/ProcessLib/BoundaryCondition/DirichletBoundaryCondition.h b/ProcessLib/BoundaryCondition/DirichletBoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..9dccf2e1820404e44f22783026e41a37ed8ee0d0
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/DirichletBoundaryCondition.h
@@ -0,0 +1,34 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESS_LIB_DIRICHLETBC_H
+#define PROCESS_LIB_DIRICHLETBC_H
+
+#include "NumLib/IndexValueVector.h"
+#include "BoundaryCondition.h"
+
+namespace ProcessLib
+{
+class DirichletBoundaryCondition : public BoundaryCondition
+{
+public:
+ virtual NumLib::IndexValueVector<GlobalIndexType> getBCValues() = 0;
+
+ void apply(const double /*t*/,
+ GlobalVector const& /*x*/,
+ GlobalMatrix& /*K*/,
+ GlobalVector& /*b*/) override final
+ {
+ // Do nothing. Dirichlet BCs are handled specially.
+ }
+};
+
+} // namespace ProcessLib
+
+#endif // PROCESS_LIB_DIRICHLETBC_H
diff --git a/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition-impl.h b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition-impl.h
new file mode 100644
index 0000000000000000000000000000000000000000..3bcf89bc256e31878849f3bb7eceb2be368e4d7b
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition-impl.h
@@ -0,0 +1,90 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "GenericNaturalBoundaryCondition.h"
+#include "ProcessLib/Utils/CreateLocalAssemblers.h"
+#include "MeshLib/MeshSearch/NodeSearch.h"
+#include "GenericNaturalBoundaryConditionLocalAssembler.h"
+
+namespace ProcessLib
+{
+template <typename BoundaryConditionData,
+ template <typename, typename, unsigned>
+ class LocalAssemblerImplementation>
+template <typename Data>
+GenericNaturalBoundaryCondition<BoundaryConditionData,
+ LocalAssemblerImplementation>::
+ GenericNaturalBoundaryCondition(
+ typename std::enable_if<
+ std::is_same<typename std::decay<BoundaryConditionData>::type,
+ typename std::decay<Data>::type>::value,
+ unsigned const>::type integration_order,
+ NumLib::LocalToGlobalIndexMap const& dof_table_bulk,
+ int const variable_id, int const component_id,
+ unsigned const global_dim,
+ std::vector<MeshLib::Element*> const& elements, Data&& data)
+ : _data(std::forward<Data>(data)), _integration_order(integration_order)
+{
+ assert(component_id <
+ static_cast<int>(dof_table_bulk.getNumberOfComponents()));
+
+ // deep copy because the neumann bc config destroys the elements.
+ std::transform(elements.begin(), elements.end(),
+ std::back_inserter(_elements),
+ std::mem_fn(&MeshLib::Element::clone));
+
+ std::vector<MeshLib::Node*> nodes = MeshLib::getUniqueNodes(_elements);
+
+ auto const& mesh_subsets =
+ dof_table_bulk.getMeshSubsets(variable_id, component_id);
+
+ // TODO extend the node intersection to all parts of mesh_subsets, i.e.
+ // to each of the MeshSubset in the mesh_subsets.
+ _mesh_subset_all_nodes.reset(
+ mesh_subsets.getMeshSubset(0).getIntersectionByNodes(nodes));
+ std::unique_ptr<MeshLib::MeshSubsets> all_mesh_subsets{
+ new MeshLib::MeshSubsets{_mesh_subset_all_nodes.get()}};
+
+ // Create local DOF table from intersected mesh subsets for the given
+ // variable and component ids.
+ _dof_table_boundary.reset(dof_table_bulk.deriveBoundaryConstrainedMap(
+ variable_id, component_id, std::move(all_mesh_subsets), _elements));
+
+ createLocalAssemblers<LocalAssemblerImplementation>(
+ global_dim, _elements, *_dof_table_boundary, _integration_order,
+ _local_assemblers, _data);
+}
+
+template <typename BoundaryConditionData,
+ template <typename, typename, unsigned>
+ class LocalAssemblerImplementation>
+GenericNaturalBoundaryCondition<
+ BoundaryConditionData,
+ LocalAssemblerImplementation>::~GenericNaturalBoundaryCondition()
+{
+ for (auto e : _elements)
+ delete e;
+}
+
+template <typename BoundaryConditionData,
+ template <typename, typename, unsigned>
+ class LocalAssemblerImplementation>
+void GenericNaturalBoundaryCondition<
+ BoundaryConditionData,
+ LocalAssemblerImplementation>::apply(const double t,
+ const GlobalVector& x,
+ GlobalMatrix& K,
+ GlobalVector& b)
+{
+ GlobalExecutor::executeMemberOnDereferenced(
+ &GenericNaturalBoundaryConditionLocalAssemblerInterface::assemble,
+ _local_assemblers, *_dof_table_boundary, t, x, K, b);
+}
+
+} // ProcessLib
diff --git a/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition.h b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..8b29a94960da5e194006e5a3dec6f40daec4c9c7
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryCondition.h
@@ -0,0 +1,77 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_GENERICNATURALBOUNDARYCONDITION_H
+#define PROCESSLIB_GENERICNATURALBOUNDARYCONDITION_H
+
+#include "MeshLib/MeshSubset.h"
+#include "BoundaryCondition.h"
+
+namespace ProcessLib
+{
+class GenericNaturalBoundaryConditionLocalAssemblerInterface;
+
+template <typename BoundaryConditionData,
+ template <typename, typename, unsigned>
+ class LocalAssemblerImplementation>
+class GenericNaturalBoundaryCondition final : public BoundaryCondition
+{
+public:
+ /// Create a boundary condition process from given config,
+ /// DOF-table, and a mesh subset for a given variable and its component.
+ /// A local DOF-table, a subset of the given one, is constructed.
+ template <typename Data>
+ GenericNaturalBoundaryCondition(
+ typename std::enable_if<
+ std::is_same<typename std::decay<BoundaryConditionData>::type,
+ typename std::decay<Data>::type>::value,
+ unsigned const>::type integration_order,
+ NumLib::LocalToGlobalIndexMap const& dof_table_bulk,
+ int const variable_id, int const component_id,
+ unsigned const global_dim,
+ std::vector<MeshLib::Element*> const& elements, Data&& data);
+
+ ~GenericNaturalBoundaryCondition() override;
+
+ /// Calls local assemblers which calculate their contributions to the global
+ /// matrix and the right-hand-side.
+ void apply(const double t,
+ GlobalVector const& x,
+ GlobalMatrix& K,
+ GlobalVector& b) override;
+
+private:
+ /// Data used in the assembly of the specific boundary condition.
+ BoundaryConditionData _data;
+
+ /// Vector of lower-dimensional elements on which the boundary condition is
+ /// defined.
+ std::vector<MeshLib::Element*> _elements;
+
+ std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_all_nodes;
+
+ /// Local dof table, a subset of the global one restricted to the
+ /// participating #_elements of the boundary condition.
+ std::unique_ptr<NumLib::LocalToGlobalIndexMap> _dof_table_boundary;
+
+ /// Integration order for integration over the lower-dimensional elements of
+ /// the #_function.
+ unsigned const _integration_order;
+
+ /// Local assemblers for each element of #_elements.
+ std::vector<
+ std::unique_ptr<GenericNaturalBoundaryConditionLocalAssemblerInterface>>
+ _local_assemblers;
+};
+
+} // ProcessLib
+
+#include "GenericNaturalBoundaryCondition-impl.h"
+
+#endif // PROCESSLIB_GENERICNATURALBOUNDARYCONDITION_H
diff --git a/ProcessLib/BoundaryCondition/GenericNaturalBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryConditionLocalAssembler.h
new file mode 100644
index 0000000000000000000000000000000000000000..0469b00f96a8870395f573cf0f5d68241be34e21
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/GenericNaturalBoundaryConditionLocalAssembler.h
@@ -0,0 +1,56 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_GENERICNATURALBOUNDARYCONDITIONLOCALASSEMBLER_H
+#define PROCESSLIB_GENERICNATURALBOUNDARYCONDITIONLOCALASSEMBLER_H
+
+#include "NumLib/Fem/ShapeMatrixPolicy.h"
+#include "ProcessLib/Utils/InitShapeMatrices.h"
+
+namespace ProcessLib
+{
+class GenericNaturalBoundaryConditionLocalAssemblerInterface
+{
+public:
+ virtual ~GenericNaturalBoundaryConditionLocalAssemblerInterface() = default;
+
+ virtual void assemble(
+ std::size_t const id,
+ NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t,
+ const GlobalVector& x, GlobalMatrix& K, GlobalVector& b) = 0;
+};
+
+template <typename ShapeFunction, typename IntegrationMethod,
+ unsigned GlobalDim>
+class GenericNaturalBoundaryConditionLocalAssembler
+ : public GenericNaturalBoundaryConditionLocalAssemblerInterface
+{
+protected:
+ using ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, GlobalDim>;
+ using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType;
+ using NodalVectorType = typename ShapeMatricesType::NodalVectorType;
+
+public:
+ GenericNaturalBoundaryConditionLocalAssembler(
+ MeshLib::Element const& e, unsigned const integration_order)
+ : _shape_matrices(initShapeMatrices<ShapeFunction, ShapeMatricesType,
+ IntegrationMethod, GlobalDim>(
+ e, integration_order)),
+ _integration_order(integration_order)
+ {
+ }
+
+protected:
+ std::vector<typename ShapeMatricesType::ShapeMatrices> _shape_matrices;
+ unsigned const _integration_order;
+};
+
+} // ProcessLib
+
+#endif // PROCESSLIB_GENERICNATURALBOUNDARYCONDITIONLOCALASSEMBLER_H
diff --git a/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.cpp b/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..fe6ce1e9bf7b86a8393972c257fd2ab5e52ecddb
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.cpp
@@ -0,0 +1,67 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "UniformDirichletBoundaryCondition.h"
+
+#include <algorithm>
+#include <vector>
+#include <logog/include/logog.hpp>
+
+#include "MeshGeoToolsLib/MeshNodeSearcher.h"
+#include "BoundaryConditionConfig.h"
+
+namespace ProcessLib
+{
+std::unique_ptr<UniformDirichletBoundaryCondition>
+createUniformDirichletBoundaryCondition(
+ BoundaryConditionConfig const& config,
+ MeshGeoToolsLib::MeshNodeSearcher& searcher,
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ int const variable_id)
+{
+ DBUG("Constructing UniformDirichletBoundaryCondition from config.");
+ //! \ogs_file_param{boundary_condition__type}
+ config.config.checkConfigParameter("type", "UniformDirichlet");
+
+ //! \ogs_file_param{boundary_condition__UniformDirichlet__value}
+ auto const value = config.config.getConfigParameter<double>("value");
+ DBUG("Using value %g", value);
+
+ // Find nodes' ids on the given mesh on which this boundary condition
+ // is defined.
+ std::vector<std::size_t> ids = searcher.getMeshNodeIDs(config.geometry);
+
+ NumLib::IndexValueVector<GlobalIndexType> bc;
+
+ // convert mesh node ids to global index for the given component
+ bc.ids.reserve(bc.ids.size() + ids.size());
+ bc.values.reserve(bc.values.size() + ids.size());
+ for (auto& id : ids) {
+ MeshLib::Location l(searcher.getMeshId(), MeshLib::MeshItemType::Node,
+ id);
+ // TODO: that might be slow, but only done once
+ const auto g_idx =
+ dof_table.getGlobalIndex(l, variable_id, config.component_id);
+ // For the DDC approach (e.g. with PETSc option), the negative
+ // index of g_idx means that the entry by that index is a ghost one,
+ // which should be dropped. Especially for PETSc routines MatZeroRows
+ // and MatZeroRowsColumns, which are called to apply the Dirichlet BC,
+ // the negative index is not accepted like other matrix or vector
+ // PETSc routines. Therefore, the following if-condition is applied.
+ if (g_idx >= 0) {
+ bc.ids.emplace_back(g_idx);
+ bc.values.emplace_back(value);
+ }
+ }
+
+ return std::unique_ptr<UniformDirichletBoundaryCondition>(
+ new UniformDirichletBoundaryCondition(std::move(bc)));
+}
+
+} // namespace ProcessLib
diff --git a/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.h b/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..f2c6e6b54a5774e3090bf1b2128ea92f9e2c787c
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/UniformDirichletBoundaryCondition.h
@@ -0,0 +1,64 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESS_LIB_BOUNDARY_CONDITION_H_
+#define PROCESS_LIB_BOUNDARY_CONDITION_H_
+
+#include "BaseLib/ConfigTree.h"
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
+#include "NumLib/NumericsConfig.h" // for GlobalIndexType
+
+#include "DirichletBoundaryCondition.h"
+
+namespace GeoLib
+{
+class GeoObject;
+}
+
+namespace MeshGeoToolsLib
+{
+class MeshNodeSearcher;
+}
+
+namespace ProcessLib
+{
+struct BoundaryConditionConfig;
+
+/// The UniformDirichletBoundaryCondition class describes a constant in space
+/// and time Dirichlet boundary condition.
+/// The expected parameter in the passed configuration is "value" which, when
+/// not present defaults to zero.
+class UniformDirichletBoundaryCondition : public DirichletBoundaryCondition
+{
+public:
+ UniformDirichletBoundaryCondition(
+ NumLib::IndexValueVector<GlobalIndexType>&& bc)
+ : _bc(std::move(bc))
+ {
+ }
+
+ NumLib::IndexValueVector<GlobalIndexType> getBCValues()
+ {
+ return std::move(_bc);
+ }
+
+private:
+ NumLib::IndexValueVector<GlobalIndexType> _bc;
+};
+
+std::unique_ptr<UniformDirichletBoundaryCondition>
+createUniformDirichletBoundaryCondition(
+ BoundaryConditionConfig const& config,
+ MeshGeoToolsLib::MeshNodeSearcher& searcher,
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ int const variable_id);
+
+} // namespace ProcessLib
+
+#endif // PROCESS_LIB_BOUNDARY_CONDITION_H_
diff --git a/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.cpp b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..d14754f29a61c4b10fe7e57482d619c9d2374cba
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.cpp
@@ -0,0 +1,41 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "UniformNeumannBoundaryCondition.h"
+#include "MeshGeoToolsLib/BoundaryElementsSearcher.h"
+#include "MeshLib/MeshSearch/NodeSearch.h"
+#include "ProcessLib/Utils/CreateLocalAssemblers.h"
+#include "BoundaryConditionConfig.h"
+
+namespace ProcessLib
+{
+std::unique_ptr<UniformNeumannBoundaryCondition>
+createUniformNeumannBoundaryCondition(
+ BoundaryConditionConfig const& config,
+ MeshGeoToolsLib::BoundaryElementsSearcher& searcher,
+ NumLib::LocalToGlobalIndexMap const& dof_table, int const variable_id,
+ unsigned const integration_order, const unsigned global_dim)
+{
+ DBUG("Constructing NeumannBcConfig from config.");
+ //! \ogs_file_param{boundary_condition__type}
+ config.config.checkConfigParameter("type", "UniformNeumann");
+
+ //! \ogs_file_param{boundary_condition__UniformNeumann__value}
+ double const value = config.config.getConfigParameter<double>("value");
+ DBUG("Using value %g", value);
+
+ auto& elems = searcher.getBoundaryElements(config.geometry);
+
+ return std::unique_ptr<UniformNeumannBoundaryCondition>(
+ new UniformNeumannBoundaryCondition(integration_order, dof_table,
+ variable_id, config.component_id,
+ global_dim, elems, value));
+}
+
+} // ProcessLib
diff --git a/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.h b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..9764f2686bc4373868cf48b9cb03c39fc24d32ad
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryCondition.h
@@ -0,0 +1,35 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITION_H
+#define PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITION_H
+
+#include "GenericNaturalBoundaryCondition.h"
+#include "UniformNeumannBoundaryConditionLocalAssembler.h"
+
+namespace MeshGeoToolsLib
+{
+class BoundaryElementsSearcher;
+}
+
+namespace ProcessLib
+{
+using UniformNeumannBoundaryCondition = GenericNaturalBoundaryCondition<
+ double, UniformNeumannBoundaryConditionLocalAssembler>;
+
+std::unique_ptr<UniformNeumannBoundaryCondition>
+createUniformNeumannBoundaryCondition(
+ BoundaryConditionConfig const& config,
+ MeshGeoToolsLib::BoundaryElementsSearcher& searcher,
+ NumLib::LocalToGlobalIndexMap const& dof_table, int const variable_id,
+ unsigned const integration_order, unsigned const global_dim);
+
+} // ProcessLib
+
+#endif // PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITION_H
diff --git a/ProcessLib/BoundaryCondition/UniformNeumannBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryConditionLocalAssembler.h
new file mode 100644
index 0000000000000000000000000000000000000000..2256c5d448567fffb3163a436f1b4ae7fae6b2cb
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/UniformNeumannBoundaryConditionLocalAssembler.h
@@ -0,0 +1,71 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITIONLOCALASSEMBLER_H
+#define PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITIONLOCALASSEMBLER_H
+
+#include "NumLib/DOF/DOFTableUtil.h"
+
+namespace ProcessLib
+{
+template <typename ShapeFunction, typename IntegrationMethod,
+ unsigned GlobalDim>
+class UniformNeumannBoundaryConditionLocalAssembler final
+ : public GenericNaturalBoundaryConditionLocalAssembler<
+ ShapeFunction, IntegrationMethod, GlobalDim>
+{
+ using Base = GenericNaturalBoundaryConditionLocalAssembler<
+ ShapeFunction, IntegrationMethod, GlobalDim>;
+
+public:
+ /// The neumann_bc_value factor is directly integrated into the local
+ /// element matrix.
+ UniformNeumannBoundaryConditionLocalAssembler(
+ MeshLib::Element const& e,
+ std::size_t const local_matrix_size,
+ unsigned const integration_order,
+ double neumann_bc_value)
+ : Base(e, integration_order),
+ _neumann_bc_value(neumann_bc_value),
+ _local_rhs(local_matrix_size)
+ {
+ }
+
+ void assemble(std::size_t const id,
+ NumLib::LocalToGlobalIndexMap const& dof_table_boundary,
+ double const t, const GlobalVector& /*x*/,
+ GlobalMatrix& /*K*/, GlobalVector& b) override
+ {
+ (void)t; // TODO time-dependent Neumann BCs
+
+ _local_rhs.setZero();
+
+ IntegrationMethod integration_method(Base::_integration_order);
+ std::size_t const n_integration_points =
+ integration_method.getNumberOfPoints();
+
+ for (std::size_t ip(0); ip < n_integration_points; ip++) {
+ auto const& sm = Base::_shape_matrices[ip];
+ auto const& wp = integration_method.getWeightedPoint(ip);
+ _local_rhs.noalias() +=
+ sm.N * _neumann_bc_value * sm.detJ * wp.getWeight();
+ }
+
+ auto const indices = NumLib::getIndices(id, dof_table_boundary);
+ b.add(indices, _local_rhs);
+ }
+
+private:
+ double const _neumann_bc_value;
+ typename Base::NodalVectorType _local_rhs;
+};
+
+} // namespace ProcessLib
+
+#endif // PROCESSLIB_UNIFORMNEUMANNBOUNDARYCONDITIONLOCALASSEMBLER_H