diff --git a/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/c_NonuniformVariableDependentNeumann.md b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/c_NonuniformVariableDependentNeumann.md
new file mode 100644
index 0000000000000000000000000000000000000000..9d198c12335b0caeba7a7b9c1fa709ae85c2dd87
--- /dev/null
+++ b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/c_NonuniformVariableDependentNeumann.md
@@ -0,0 +1,3 @@
+Declares a Neumann boundary condition for a two variable process that is non-uniform in space and depending on both variable values on the mesh nodes of the boundary mesh.
+The flux consists of 4 summands: 1) a constant flux, 2) a flux depending on the process variables value at the mesh node, 3) a flux depending on the other variables value at the mesh node and 3) a flux depending on the product of both variables at the mesh node.
+If the Newton-Raphson method is applied for the global system, please select the numerical Jacobian calculation.
diff --git a/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_current_variable_name.md b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_current_variable_name.md
new file mode 100644
index 0000000000000000000000000000000000000000..1bdf42c51d49c02df3c99f4047873e108daad705
--- /dev/null
+++ b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_current_variable_name.md
@@ -0,0 +1,3 @@
+This specifies a nodal field defined on the surface mesh.
+
+From the nodal field the coefficient values, which are multiplied by the process variable to calculate the current variable dependent flux, which is applied via this Neumann BC, are read.
diff --git a/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_mixed_variables_name.md b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_mixed_variables_name.md
new file mode 100644
index 0000000000000000000000000000000000000000..828ed39873b07c1d2c3937b6adcf80f53a9171fa
--- /dev/null
+++ b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_mixed_variables_name.md
@@ -0,0 +1,3 @@
+This specifies a nodal field defined on the surface mesh.
+
+From the nodal field the coefficient values, which are multiplied by the product of both variables to calculate the mixed variable dependent flux, which is applied via this Neumann BC, are read.
diff --git a/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_other_variable_name.md b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_other_variable_name.md
new file mode 100644
index 0000000000000000000000000000000000000000..bb42d4abbe625c5cdcdabf6228c3f853cd8c833e
--- /dev/null
+++ b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_coefficient_other_variable_name.md
@@ -0,0 +1,3 @@
+This specifies a nodal field defined on the surface mesh.
+
+From the nodal field the coefficient values, which are multiplied by the other variable to calculate the other variable dependent flux, which is applied via this Neumann BC, are read.
diff --git a/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_constant_name.md b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_constant_name.md
new file mode 100644
index 0000000000000000000000000000000000000000..24c508339619f7a4eb2c2d9c8ba08d9fe1d4092a
--- /dev/null
+++ b/Documentation/ProjectFile/prj/process_variables/process_variable/boundary_conditions/boundary_condition/NonuniformVariableDependentNeumann/t_constant_name.md
@@ -0,0 +1,4 @@
+This specifies a nodal field defined on the surface mesh.
+
+From the nodal field the constant flux values, which are applied via this Neumann BC,
+are read.
diff --git a/MeshLib/Properties-impl.h b/MeshLib/Properties-impl.h
index 6af4d52bd12384bcf5019beab0b5862ca7ba8d21..960da18e007bbee9cd9682e8a08671ee2c60e07f 100644
--- a/MeshLib/Properties-impl.h
+++ b/MeshLib/Properties-impl.h
@@ -20,7 +20,7 @@ PropertyVector<T>* Properties::createNewPropertyVector(
_properties.find(name)
);
if (it != _properties.end()) {
- ERR("A property of the name \"%s\" is already assigned to the mesh.",
+ ERR("A property of the name '%s' is already assigned to the mesh.",
name.c_str());
return nullptr;
}
@@ -48,7 +48,7 @@ PropertyVector<T>* Properties::createNewPropertyVector(
_properties.find(name)
);
if (it != _properties.end()) {
- ERR("A property of the name \"%s\" already assigned to the mesh.",
+ ERR("A property of the name '%s' already assigned to the mesh.",
name.c_str());
return nullptr;
}
@@ -128,3 +128,33 @@ PropertyVector<T>* Properties::getPropertyVector(std::string const& name)
return dynamic_cast<PropertyVector<T>*>(it->second);
}
+template <typename T>
+PropertyVector<T> const* Properties::getPropertyVector(
+ std::string const& name, MeshItemType const item_type,
+ int const n_components) const
+{
+ auto const it = _properties.find(name);
+ if (it == _properties.end())
+ {
+ OGS_FATAL("A property with name '%s' does not exist.", name.c_str());
+ }
+
+ auto property = dynamic_cast<PropertyVector<T>*>(it->second);
+ if (property == nullptr)
+ {
+ OGS_FATAL("Could not cast property '%s' to given type.", name.c_str());
+ }
+ if (property->getMeshItemType() != item_type)
+ {
+ OGS_FATAL(
+ "The PropertyVector '%s' has a different type than requested. A "
+ "'%s' field is requested.",
+ name.c_str(), toString(item_type));
+ }
+ if (property->getNumberOfComponents() != n_components)
+ {
+ OGS_FATAL("'%s' does not have the right number of components.",
+ name.c_str());
+ }
+ return property;
+}
diff --git a/MeshLib/Properties.h b/MeshLib/Properties.h
index 1ba2a4f53231eadd2aa42d7bcf5cbb5db4edb6db..220889273ce304535c789317a54ffc5c6ab45ab5 100644
--- a/MeshLib/Properties.h
+++ b/MeshLib/Properties.h
@@ -95,6 +95,14 @@ public:
template <typename T>
PropertyVector<T>* getPropertyVector(std::string const& name);
+ /// Returns a property vector with given \c name, \c item_type and \c
+ /// number_of_components or aborts calling OGS_FATAL if no such property
+ /// vector exists.
+ template <typename T>
+ PropertyVector<T> const* getPropertyVector(std::string const& name,
+ MeshItemType const item_type,
+ int const n_components) const;
+
void removePropertyVector(std::string const& name);
/// Check if a PropertyVector accessible by the name is already
diff --git a/ProcessLib/BoundaryCondition/CreateBoundaryCondition.cpp b/ProcessLib/BoundaryCondition/CreateBoundaryCondition.cpp
index 3bc07209d7211647e88615c3ad6a390b36ad56b0..98b5a99c1e6a8e863f2c488efd17d79fc3891b01 100644
--- a/ProcessLib/BoundaryCondition/CreateBoundaryCondition.cpp
+++ b/ProcessLib/BoundaryCondition/CreateBoundaryCondition.cpp
@@ -16,6 +16,7 @@
#include "NeumannBoundaryCondition.h"
#include "NonuniformDirichletBoundaryCondition.h"
#include "NonuniformNeumannBoundaryCondition.h"
+#include "NonuniformVariableDependentNeumannBoundaryCondition.h"
#include "NormalTractionBoundaryCondition.h"
#include "PhaseFieldIrreversibleDamageOracleBoundaryCondition.h"
#include "RobinBoundaryCondition.h"
@@ -80,6 +81,16 @@ std::unique_ptr<BoundaryCondition> createBoundaryCondition(
*config.component_id, integration_order, shapefunction_order,
bulk_mesh);
}
+
+ if (type == "NonuniformVariableDependentNeumann")
+ {
+ return ProcessLib::
+ createNonuniformVariableDependentNeumannBoundaryCondition(
+ config.config, config.boundary_mesh, dof_table, variable_id,
+ *config.component_id, integration_order, shapefunction_order,
+ bulk_mesh);
+ }
+
if (type == "Python")
{
#ifdef OGS_USE_PYTHON
diff --git a/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.cpp b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..5e6829484fb859fc49690c95de82ee8ba5e65ce0
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.cpp
@@ -0,0 +1,101 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2018, 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 "NonuniformVariableDependentNeumannBoundaryCondition.h"
+
+#include "MeshLib/IO/readMeshFromFile.h"
+#include "ProcessLib/Utils/ProcessUtils.h"
+
+namespace ProcessLib
+{
+std::unique_ptr<NonuniformVariableDependentNeumannBoundaryCondition>
+createNonuniformVariableDependentNeumannBoundaryCondition(
+ BaseLib::ConfigTree const& config, MeshLib::Mesh const& boundary_mesh,
+ NumLib::LocalToGlobalIndexMap const& dof_table, int const variable_id,
+ int const component_id, unsigned const integration_order,
+ unsigned const shapefunction_order, MeshLib::Mesh const& bulk_mesh)
+{
+ DBUG("Constructing NonuniformVariableDependentNeumann BC from config.");
+ //! \ogs_file_param{prj__process_variables__process_variable__boundary_conditions__boundary_condition__type}
+ config.checkConfigParameter("type", "NonuniformVariableDependentNeumann");
+ if (dof_table.getNumberOfVariables() != 2)
+ {
+ OGS_FATAL(
+ "NonuniformVariableDependentNeumann BC only implemented for 2 "
+ "variable processes.");
+ }
+ assert(variable_id == 0 || variable_id == 1);
+
+ // TODO finally use ProcessLib::Parameter here
+ auto const constant_name =
+ //! \ogs_file_param{prj__process_variables__process_variable__boundary_conditions__boundary_condition__NonuniformVariableDependentNeumann__constant_name}
+ config.getConfigParameter<std::string>("constant_name");
+ auto const& constant =
+ *boundary_mesh.getProperties().getPropertyVector<double>(
+ constant_name, MeshLib::MeshItemType::Node, 1);
+
+ auto const coefficient_current_variable_name =
+ //! \ogs_file_param{prj__process_variables__process_variable__boundary_conditions__boundary_condition__NonuniformVariableDependentNeumann__coefficient_current_variable_name}
+ config.getConfigParameter<std::string>(
+ "coefficient_current_variable_name");
+ auto const& coefficient_current_variable =
+ *boundary_mesh.getProperties().getPropertyVector<double>(
+ coefficient_current_variable_name, MeshLib::MeshItemType::Node, 1);
+
+ auto const coefficient_other_variable_name =
+ //! \ogs_file_param{prj__process_variables__process_variable__boundary_conditions__boundary_condition__NonuniformVariableDependentNeumann__coefficient_other_variable_name}
+ config.getConfigParameter<std::string>(
+ "coefficient_other_variable_name");
+ auto const& coefficient_other_variable =
+ *boundary_mesh.getProperties().getPropertyVector<double>(
+ coefficient_other_variable_name, MeshLib::MeshItemType::Node, 1);
+
+ auto const coefficient_mixed_variables_name =
+ //! \ogs_file_param{prj__process_variables__process_variable__boundary_conditions__boundary_condition__NonuniformVariableDependentNeumann__coefficient_mixed_variables_name}
+ config.getConfigParameter<std::string>(
+ "coefficient_mixed_variables_name");
+ auto const& coefficient_mixed_variables =
+ *boundary_mesh.getProperties().getPropertyVector<double>(
+ coefficient_mixed_variables_name, MeshLib::MeshItemType::Node, 1);
+
+ std::string const mapping_to_bulk_nodes_property = "bulk_node_ids";
+ boundary_mesh.getProperties().getPropertyVector<std::size_t>(
+ mapping_to_bulk_nodes_property, MeshLib::MeshItemType::Node, 1);
+
+ std::vector<MeshLib::Node*> const& bc_nodes = boundary_mesh.getNodes();
+ MeshLib::MeshSubset bc_mesh_subset(boundary_mesh, bc_nodes);
+ auto const& dof_table_boundary_other_variable =
+ *dof_table.deriveBoundaryConstrainedMap(
+ (variable_id + 1) % 2, {component_id}, std::move(bc_mesh_subset));
+
+ // In case of partitioned mesh the boundary could be empty, i.e. there is no
+ // boundary condition.
+#ifdef USE_PETSC
+ // This can be extracted to createBoundaryCondition() but then the config
+ // parameters are not read and will cause an error.
+ // TODO (naumov): Add a function to ConfigTree for skipping the tags of the
+ // subtree and move the code up in createBoundaryCondition().
+ if (boundary_mesh.getDimension() == 0 &&
+ boundary_mesh.getNumberOfNodes() == 0 &&
+ boundary_mesh.getNumberOfElements() == 0)
+ {
+ return nullptr;
+ }
+#endif // USE_PETSC
+
+ return std::make_unique<
+ NonuniformVariableDependentNeumannBoundaryCondition>(
+ integration_order, shapefunction_order, dof_table, variable_id,
+ component_id, bulk_mesh.getDimension(), boundary_mesh,
+ NonuniformVariableDependentNeumannBoundaryConditionData{
+ constant, coefficient_current_variable, coefficient_other_variable,
+ coefficient_mixed_variables, dof_table_boundary_other_variable});
+}
+
+} // namespace ProcessLib
diff --git a/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.h b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.h
new file mode 100644
index 0000000000000000000000000000000000000000..a9b8b57885e1c27648640cacc215a293669638a8
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryCondition.h
@@ -0,0 +1,30 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include "GenericNonuniformNaturalBoundaryCondition.h"
+#include "MeshLib/PropertyVector.h"
+#include "NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler.h"
+
+namespace ProcessLib
+{
+using NonuniformVariableDependentNeumannBoundaryCondition =
+ GenericNonuniformNaturalBoundaryCondition<
+ NonuniformVariableDependentNeumannBoundaryConditionData,
+ NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler>;
+
+std::unique_ptr<NonuniformVariableDependentNeumannBoundaryCondition>
+createNonuniformVariableDependentNeumannBoundaryCondition(
+ BaseLib::ConfigTree const& config, MeshLib::Mesh const& boundary_mesh,
+ NumLib::LocalToGlobalIndexMap const& dof_table, int const variable_id,
+ int const component_id, unsigned const integration_order,
+ unsigned const shapefunction_order, const MeshLib::Mesh& bulk_mesh);
+
+} // namespace ProcessLib
diff --git a/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler.h
new file mode 100644
index 0000000000000000000000000000000000000000..4bcae69fa269b156730c9f6e328a0c8fd60e91ce
--- /dev/null
+++ b/ProcessLib/BoundaryCondition/NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler.h
@@ -0,0 +1,128 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include "MeshLib/PropertyVector.h"
+#include "NumLib/DOF/DOFTableUtil.h"
+#include "NumLib/Function/Interpolation.h"
+#include "ProcessLib/Parameter/MeshNodeParameter.h"
+
+#include "GenericNonuniformNaturalBoundaryConditionLocalAssembler.h"
+
+namespace ProcessLib
+{
+struct NonuniformVariableDependentNeumannBoundaryConditionData
+{
+ /* TODO use Parameter */
+ MeshLib::PropertyVector<double> const& constant;
+ MeshLib::PropertyVector<double> const& coefficient_current_variable;
+ MeshLib::PropertyVector<double> const& coefficient_other_variable;
+ MeshLib::PropertyVector<double> const& coefficient_mixed_variables;
+ // Used for mapping boundary nodes to bulk nodes.
+ NumLib::LocalToGlobalIndexMap const& dof_table_boundary_other_variable;
+};
+
+template <typename ShapeFunction, typename IntegrationMethod,
+ unsigned GlobalDim>
+class NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler final
+ : public GenericNonuniformNaturalBoundaryConditionLocalAssembler<
+ ShapeFunction, IntegrationMethod, GlobalDim>
+{
+ using Base = GenericNonuniformNaturalBoundaryConditionLocalAssembler<
+ ShapeFunction, IntegrationMethod, GlobalDim>;
+ using NodalVectorType = typename Base::NodalVectorType;
+ using NodalMatrixType = typename Base::NodalMatrixType;
+
+public:
+ /// The neumann_bc_term factor is directly integrated into the local
+ /// element matrix.
+ NonuniformVariableDependentNeumannBoundaryConditionLocalAssembler(
+ MeshLib::Element const& e,
+ std::size_t const local_matrix_size,
+ bool const is_axially_symmetric,
+ unsigned const integration_order,
+ NonuniformVariableDependentNeumannBoundaryConditionData const& data)
+ : Base(e, is_axially_symmetric, integration_order),
+ _data(data),
+ _local_matrix_size(local_matrix_size)
+ {
+ }
+
+ void assemble(std::size_t const mesh_item_id,
+ NumLib::LocalToGlobalIndexMap const& dof_table_boundary,
+ double const t, const GlobalVector& x, GlobalMatrix& /*K*/,
+ GlobalVector& b, GlobalMatrix* /*Jac*/) override
+ {
+ NodalVectorType _local_rhs(_local_matrix_size);
+ _local_rhs.setZero();
+ MeshNodeParameter<double> constant_values{"ConstantValues",
+ _data.constant};
+ MeshNodeParameter<double> coefficient_current_variable_values{
+ "CurrentVariableValues", _data.coefficient_current_variable};
+ MeshNodeParameter<double> coefficient_other_variable_values{
+ "OtherVariableValues", _data.coefficient_other_variable};
+ MeshNodeParameter<double> coefficient_mixed_variables_values{
+ "MixedVariablesValues", _data.coefficient_mixed_variables};
+ // Get element nodes for the interpolation from nodes to
+ // integration point.
+ NodalVectorType const constant_node_values =
+ constant_values.getNodalValuesOnElement(Base::_element, t);
+ NodalVectorType const coefficient_current_variable_node_values =
+ coefficient_current_variable_values.getNodalValuesOnElement(
+ Base::_element, t);
+ NodalVectorType const coefficient_other_variable_node_values =
+ coefficient_other_variable_values.getNodalValuesOnElement(
+ Base::_element, t);
+ NodalVectorType const coefficient_mixed_variables_node_values =
+ coefficient_mixed_variables_values.getNodalValuesOnElement(
+ Base::_element, t);
+ unsigned const n_integration_points = _local_matrix_size;
+
+ auto const indices_current_variable =
+ NumLib::getIndices(mesh_item_id, dof_table_boundary);
+ auto const indices_other_variable = NumLib::getIndices(
+ mesh_item_id, _data.dof_table_boundary_other_variable);
+ std::vector<double> const local_current_variable =
+ x.get(indices_current_variable);
+ std::vector<double> const local_other_variable =
+ x.get(indices_other_variable);
+
+ for (unsigned ip = 0; ip < n_integration_points; ip++)
+ {
+ auto const& n_and_weight = Base::_ns_and_weights[ip];
+ auto const& N = n_and_weight.N;
+ auto const& w = n_and_weight.weight;
+
+ double current_variable_int_pt = 0.0;
+ double other_variable_int_pt = 0.0;
+
+ NumLib::shapeFunctionInterpolate(local_current_variable, N,
+ current_variable_int_pt);
+ NumLib::shapeFunctionInterpolate(local_other_variable, N,
+ other_variable_int_pt);
+ NodalVectorType const neumann_node_values =
+ constant_node_values +
+ coefficient_current_variable_node_values *
+ current_variable_int_pt +
+ coefficient_other_variable_node_values * other_variable_int_pt +
+ coefficient_mixed_variables_node_values *
+ current_variable_int_pt * other_variable_int_pt;
+ _local_rhs.noalias() += N * neumann_node_values.dot(N) * w;
+ }
+
+ b.add(indices_current_variable, _local_rhs);
+ }
+
+private:
+ NonuniformVariableDependentNeumannBoundaryConditionData const& _data;
+ std::size_t const _local_matrix_size;
+};
+
+} // namespace ProcessLib
diff --git a/ProcessLib/ComponentTransport/Tests.cmake b/ProcessLib/ComponentTransport/Tests.cmake
index d205bafce1d2e111b1fde8cbf5a2ce9c9f8ee69e..c0075687e6764c9ae877d6894649ad245d4c6f4e 100644
--- a/ProcessLib/ComponentTransport/Tests.cmake
+++ b/ProcessLib/ComponentTransport/Tests.cmake
@@ -395,3 +395,29 @@ AddTest(
3D1P-GWFlow_1_reference.vtu out_ogs5_H_pcs_0_ts_1_t_10000000.000000.vtu NODAL_VELOCITY1 darcy_velocity 1e-10 1.4e-2
VIS out_ogs5_H_pcs_0_ts_1_t_10000000.000000.vtu
)
+
+
+AddTest(
+ NAME 1D_ComponentTransport_VariableDependentBoundary
+ PATH Parabolic/ComponentTransport/VariableNeumannBoundary
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS vdbc_input.prj
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ DIFF_DATA
+ vdbc_pcs_0_ts_0_t_0.000000_expected.vtu vdbc_pcs_0_ts_0_t_0.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu vdbc_pcs_0_ts_1590_t_6000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu vdbc_pcs_0_ts_3990_t_30000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu vdbc_pcs_0_ts_9990_t_90000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu vdbc_pcs_0_ts_15990_t_150000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu vdbc_pcs_0_ts_21990_t_210000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu vdbc_pcs_0_ts_25990_t_250000.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_0_t_0.000000_expected.vtu vdbc_pcs_0_ts_0_t_0.000000.vtu pressure pressure 1e-5 1e-4
+ vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu vdbc_pcs_0_ts_1590_t_6000.000000.vtu concentration concentration 1e-5 1e-4
+ vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu vdbc_pcs_0_ts_3990_t_30000.000000.vtu concentration concentration 1e-5 1e-4
+ vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu vdbc_pcs_0_ts_9990_t_90000.000000.vtu concentration concentration 1e-5 1e-4
+ vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu vdbc_pcs_0_ts_15990_t_150000.000000.vtu concentration concentration 1e-5 1e-4
+ vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu vdbc_pcs_0_ts_21990_t_210000.000000.vtu concentration concentration 1e-5 1e-4
+ vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu vdbc_pcs_0_ts_25990_t_250000.000000.vtu concentration concentration 1e-5 1e-4
+)
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.prj b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.prj
new file mode 100644
index 0000000000000000000000000000000000000000..27a339d9fc366442bfad8312cbbc1b4ac2f126d1
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.prj
@@ -0,0 +1,238 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProject>
+ <meshes>
+ <mesh>vdbc_input.vtu</mesh>
+ <mesh>vdbc_input_rightBoundary.vtu</mesh>
+ <mesh>vdbc_input_leftBoundary.vtu</mesh>
+ </meshes>
+ <processes>
+ <process>
+ <name>hc</name>
+ <type>ComponentTransport</type>
+ <integration_order>2</integration_order>
+ <process_variables>
+ <concentration>concentration</concentration>
+ <pressure>pressure</pressure>
+ </process_variables>
+ <fluid>
+ <density>
+ <type>Constant</type>
+ <value>1000</value>
+ </density>
+ <viscosity>
+ <type>Constant</type>
+ <value>1.0e-3</value>
+ </viscosity>
+ </fluid>
+ <porous_medium>
+ <porous_medium id="0">
+ <permeability>
+ <permeability_tensor_entries>kappa1</permeability_tensor_entries>
+ <type>Constant</type>
+ </permeability>
+ <porosity>
+ <type>Constant</type>
+ <porosity_parameter>constant_porosity_parameter</porosity_parameter>
+ </porosity>
+ <storage>
+ <type>Constant</type>
+ <value>0</value>
+ </storage>
+ </porous_medium>
+ </porous_medium>
+ <fluid_reference_density>rho_fluid</fluid_reference_density>
+ <solute_dispersivity_longitudinal>beta_l</solute_dispersivity_longitudinal>
+ <solute_dispersivity_transverse>beta_t</solute_dispersivity_transverse>
+ <molecular_diffusion_coefficient>Dm</molecular_diffusion_coefficient>
+ <retardation_factor>retardation</retardation_factor>
+ <decay_rate>decay</decay_rate>
+ <specific_body_force>0 0</specific_body_force>
+ <secondary_variables>
+ <secondary_variable type="static" internal_name="darcy_velocity" output_name="darcy_velocity"/>
+ </secondary_variables>
+ </process>
+ </processes>
+ <time_loop>
+ <processes>
+ <process ref="hc">
+ <nonlinear_solver>basic_picard</nonlinear_solver>
+ <convergence_criterion>
+ <type>PerComponentDeltaX</type>
+ <norm_type>NORM2</norm_type>
+ <reltols>5e-8 5e-8</reltols>
+ </convergence_criterion>
+ <time_discretization>
+ <type>BackwardEuler</type>
+ </time_discretization>
+ <output>
+ <variables>
+ <variable>concentration</variable>
+ <variable>pressure</variable>
+ <variable>darcy_velocity</variable>
+ </variables>
+ </output>
+ <time_stepping>
+ <type>FixedTimeStepping</type>
+ <t_initial>0.0</t_initial>
+ <t_end>250000.0</t_end>
+ <timesteps>
+ <pair>
+ <repeat>500</repeat>
+ <delta_t>0.1</delta_t>
+ </pair>
+ <pair>
+ <repeat>550</repeat>
+ <delta_t>1.0</delta_t>
+ </pair>
+ <pair>
+ <repeat>540</repeat>
+ <delta_t>10.0</delta_t>
+ </pair>
+ <pair>
+ <repeat>42900</repeat>
+ <delta_t>10.0</delta_t>
+ </pair>
+ </timesteps>
+ </time_stepping>
+ </process>
+ </processes>
+ <output>
+ <type>VTK</type>
+ <prefix>vdbc</prefix>
+ <timesteps>
+ <pair>
+ <repeat>1</repeat>
+ <each_steps>1590</each_steps>
+ </pair>
+ <pair>
+ <repeat>1</repeat>
+ <each_steps>2400</each_steps>
+ </pair>
+ <pair>
+ <repeat>7000</repeat>
+ <each_steps>6000</each_steps>
+ </pair>
+ </timesteps>
+ </output>
+ </time_loop>
+ <parameters>
+ <parameter>
+ <name>rho_fluid</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </parameter>
+ <parameter>
+ <name>Dm</name>
+ <type>Constant</type>
+ <value>0</value>
+ </parameter>
+ <parameter>
+ <name>retardation</name>
+ <type>Constant</type>
+ <value>1</value>
+ </parameter>
+ <parameter>
+ <name>decay</name>
+ <type>Constant</type>
+ <value>0</value>
+ </parameter>
+ <parameter>
+ <name>beta_l</name>
+ <type>Constant</type>
+ <value>0</value>
+ </parameter>
+ <parameter>
+ <name>beta_t</name>
+ <type>Constant</type>
+ <value>0</value>
+ </parameter>
+ <parameter>
+ <name>c_ini</name>
+ <type>MeshNode</type>
+ <field_name>c_ini</field_name>
+ </parameter>
+ <parameter>
+ <name>p_ini</name>
+ <type>MeshNode</type>
+ <field_name>p_ini</field_name>
+ </parameter>
+ <parameter>
+ <name>constant_porosity_parameter</name>
+ <type>Constant</type>
+ <value>0.1</value>
+ </parameter>
+ <parameter>
+ <name>kappa1</name>
+ <type>Constant</type>
+ <values>1.e-9 0 0 1.e-9</values>
+ </parameter>
+ </parameters>
+ <process_variables>
+ <process_variable>
+ <name>concentration</name>
+ <components>1</components>
+ <order>1</order>
+ <initial_condition>c_ini</initial_condition>
+ <boundary_conditions>
+ <boundary_condition>
+ <type>NonuniformVariableDependentNeumann</type>
+ <mesh>vdbc_input_rightBoundary</mesh>
+ <constant_name>Cconstant</constant_name>
+ <coefficient_current_variable_name>Cprefac1</coefficient_current_variable_name>
+ <coefficient_other_variable_name>Cprefac2</coefficient_other_variable_name>
+ <coefficient_mixed_variables_name>Cprefac3</coefficient_mixed_variables_name>
+ </boundary_condition>
+ <boundary_condition>
+ <type>NonuniformDirichlet</type>
+ <mesh>vdbc_input_leftBoundary</mesh>
+ <field_name>c_ini</field_name>
+ </boundary_condition>
+ </boundary_conditions>
+ </process_variable>
+ <process_variable>
+ <name>pressure</name>
+ <components>1</components>
+ <order>1</order>
+ <initial_condition>p_ini</initial_condition>
+ <boundary_conditions>
+ <boundary_condition>
+ <type>NonuniformVariableDependentNeumann</type>
+ <mesh>vdbc_input_rightBoundary</mesh>
+ <constant_name>Pconstant</constant_name>
+ <coefficient_current_variable_name>Pprefac1</coefficient_current_variable_name>
+ <coefficient_other_variable_name>Pprefac2</coefficient_other_variable_name>
+ <coefficient_mixed_variables_name>Pprefac3</coefficient_mixed_variables_name>
+ </boundary_condition>
+ <boundary_condition>
+ <type>NonuniformDirichlet</type>
+ <mesh>vdbc_input_leftBoundary</mesh>
+ <field_name>p_ini</field_name>
+ </boundary_condition>
+ </boundary_conditions>
+ </process_variable>
+ </process_variables>
+ <nonlinear_solvers>
+ <nonlinear_solver>
+ <name>basic_picard</name>
+ <type>Picard</type>
+ <max_iter>200</max_iter>
+ <linear_solver>general_linear_solver</linear_solver>
+ </nonlinear_solver>
+ </nonlinear_solvers>
+ <linear_solvers>
+ <linear_solver>
+ <name>general_linear_solver</name>
+ <lis>-i bicgstab -p ilut -tol 1e-12 -maxiter 200000</lis>
+ <eigen>
+ <solver_type>BiCGSTAB</solver_type>
+ <precon_type>DIAGONAL</precon_type>
+ <max_iteration_step>10000</max_iteration_step>
+ <error_tolerance>1e-12</error_tolerance>
+ </eigen>
+ <petsc>
+ <prefix>hc</prefix>
+ <parameters>-hc_ksp_type bcgs -hc_pc_type bjacobi -hc_ksp_rtol 1e-10 -hc_ksp_max_it 20000</parameters>
+ </petsc>
+ </linear_solver>
+ </linear_solvers>
+</OpenGeoSysProject>
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..49ca50f19c9e7d9e89078a1f7365a25f67bbb038
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:af547c3c28c0bb522cd0973b0ff1b3209cd49f5ead6982ec1bdd6ceff6c0fb2e
+size 2341
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_leftBoundary.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_leftBoundary.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..bd6324ea80f2f0c2589d418081d092f08878967e
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_leftBoundary.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:9f0eee7d079e1198c49ed81f2b8a298f4b810bbfe446ff7158296ae6bf273b51
+size 3844
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_rightBoundary.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_rightBoundary.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..deb34e0468bd87d5902b8b1e834849cb87881e17
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_input_rightBoundary.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:da12f42866c85e8f1e5d008996bc18521ee0cd1e7b0458668b64a897c1ff30d0
+size 3844
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_expected.pvd b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_expected.pvd
new file mode 100644
index 0000000000000000000000000000000000000000..0bcbc8a0a039852394142fcbbff0801b66e57f87
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_expected.pvd
@@ -0,0 +1,12 @@
+<?xml version="1.0"?>
+<VTKFile type="Collection" version="0.1" byte_order="LittleEndian" compressor="vtkZLibDataCompressor">
+ <Collection>
+ <DataSet timestep="0" group="" part="0" file="vdbc_pcs_0_ts_0_t_0.000000_expected.vtu"/>
+ <DataSet timestep="6000" group="" part="0" file="vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu"/>
+ <DataSet timestep="30000" group="" part="0" file="vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu"/>
+ <DataSet timestep="90000" group="" part="0" file="vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu"/>
+ <DataSet timestep="150000" group="" part="0" file="vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu"/>
+ <DataSet timestep="210000" group="" part="0" file="vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu"/>
+ <DataSet timestep="250000" group="" part="0" file="vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu"/>
+ </Collection>
+</VTKFile>
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_0_t_0.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_0_t_0.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..1db17b1351c74829831a083c3f9500e628da59c9
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_0_t_0.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:8dd72e2d28c8f339fbf9fbd58b8199c5474ae9aa96c5fd492dfad6157831f817
+size 3244
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..ad5bd5cbd95e0012657b84fb05a1c70109945f68
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_1590_t_6000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:ee350ce13b70c94c136203119fd4a499e2338b59de2f031a0c3ea36ade59c917
+size 3521
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..3633ee484e78e0fc01c498453fab713bf76ff7cb
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:be86273dc6a463d82535edb10772fb580376f53e0eed47824464aa9eedf96ba9
+size 3548
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..d966c768619a861f7258703d6fd6cdfd6353c742
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:3e50caa3c7741c0737acf5d7d3174b36fcd5a1a8fc8c051ecaa2b38cb7d0900f
+size 3545
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..d0f2a4ec4d96222b3d60eff457962d7dabbd0395
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:0a58d8b939bd0a48355e139a8e23da4a1b96c4b8a4871394df445cf90005c063
+size 3558
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..de8e54d515e3307ab1cabc72d92793c0e1d05bab
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_3990_t_30000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:ace54aff779664d34741dccfa344e5bb8bf50b5fa20babe72dc536bce4a98a58
+size 3540
diff --git a/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..dd579672c33de06eac43d8cf392b7559edc262e7
--- /dev/null
+++ b/Tests/Data/Parabolic/ComponentTransport/VariableNeumannBoundary/vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:f1e30f935c0bafa24ff7c15f0604e15a527a3946c7b0a2ce1be9df55ec0719e6
+size 3548
diff --git a/web/content/docs/benchmarks/hydro-component/HC-VDBCTest.pdf b/web/content/docs/benchmarks/hydro-component/HC-VDBCTest.pdf
new file mode 100644
index 0000000000000000000000000000000000000000..a4c1c7371d0b94fe954b342a4bd145443e0e4c3f
--- /dev/null
+++ b/web/content/docs/benchmarks/hydro-component/HC-VDBCTest.pdf
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:8a2637ce2c101de5ce27864a8d8322bd648b110ec8101c948fdf1f0ab65b4d92
+size 204834
diff --git a/web/content/docs/benchmarks/hydro-component/VDBC_num_ana_comp.png b/web/content/docs/benchmarks/hydro-component/VDBC_num_ana_comp.png
new file mode 100644
index 0000000000000000000000000000000000000000..6d23d0fa9a27440d63bf24c9a8624543af36191f
--- /dev/null
+++ b/web/content/docs/benchmarks/hydro-component/VDBC_num_ana_comp.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:e5d6e9b0ddb8d632c950e48f3e766961fae21cfba7c788678a46eb5db0caaff5
+size 84999
diff --git a/web/content/docs/benchmarks/hydro-component/vdbc.pandoc b/web/content/docs/benchmarks/hydro-component/vdbc.pandoc
new file mode 100644
index 0000000000000000000000000000000000000000..57cc9a5b92a5e101c2848ba258acdcc94d20415e
--- /dev/null
+++ b/web/content/docs/benchmarks/hydro-component/vdbc.pandoc
@@ -0,0 +1,30 @@
++++
+author = "Jasper Bathmann, Dimitri Yu. Naumov, Marc Walther"
+weight = 142
+project = "Parabolic/ComponentTransport/VariableNeumannBoundary/"
+date = "2018-11-01T10:41:09+01:00"
+title = "Variable Dependent Boundary Condition"
+
+[menu]
+ [menu.benchmarks]
+ parent = "hydro-component"
+
++++
+
+{{< data-link >}}
+
+
+## Overview
+
+The component transport process is used for the benchmark setup. Here, a analytical solution of a simple setup is derived and compared to the numerical results.
+This Benchmark is described in [this PDF](../HC-VDBCTest.pdf).
+
+For the setup and parameterization, see the chapter "Density dependent flow - The Goswami Problem" in Kolditz et al. (2012).
+
+
+## Results
+
+{{< img src="../VDBC_num_ana_comp.png" title="UPPER PART: Analytical solution on the right boundary in dependance of time $t$ of the problem indicated with red dashed line in comparison to numerical solution indicated by blue crosses; LOWER PART: development of relative error in dependance of time $t$. Grid spacing for simulations: 0.1; widest timestep 10. The relative error is below $5 \times 10^{-5}$ for all simulation times.">}}
+
+
+