Merge pull request #2061 from TomFischer/NodalSourceTermsModification

[PL] Use parameter for defining source terms.
ufz · Feb 1, 2018 · 5c74834 · 5c74834
2 parents 060810d + b36ac7d
commit 5c74834
Show file tree

Hide file tree

Showing 21 changed files with 182 additions and 47 deletions.
diff --git a/...rocess_variables/process_variable/source_terms/source_term/Nodal/t_parameter.md b/...rocess_variables/process_variable/source_terms/source_term/Nodal/t_parameter.md
@@ -0,0 +1,2 @@
+The name of the parameter that defines value that should be used for the source
+term.
diff --git a/...rj/process_variables/process_variable/source_terms/source_term/Nodal/t_value.md b/...rj/process_variables/process_variable/source_terms/source_term/Nodal/t_value.md
diff --git a/ProcessLib/Process.cpp b/ProcessLib/Process.cpp
@@ -47,6 +47,16 @@ Process::Process(
               pcs_BCs.emplace_back(BoundaryConditionCollection(parameters));
           }
           return pcs_BCs;
+      }(_process_variables.size())),
+      _source_term_collections([&](const std::size_t number_of_processes)
+                               -> std::vector<SourceTermCollection> {
+          std::vector<SourceTermCollection> pcs_sts;
+          pcs_sts.reserve(number_of_processes);
+          for (std::size_t i = 0; i < number_of_processes; i++)
+          {
+              pcs_sts.emplace_back(SourceTermCollection(parameters));
+          }
+          return pcs_sts;
       }(_process_variables.size()))
 {
 }
@@ -60,15 +70,9 @@ void Process::initializeProcessBoundaryConditionsAndSourceTerms(
     per_process_BCs.addBCsForProcessVariables(per_process_variables, dof_table,
                                               _integration_order);
 
-    std::vector<std::unique_ptr<NodalSourceTerm>> per_process_source_terms;
-    for (auto& pv : per_process_variables)
-    {
-        auto sts = pv.get().createSourceTerms(dof_table, 0, _integration_order);
-
-        std::move(sts.begin(), sts.end(),
-                  std::back_inserter(per_process_source_terms));
-    }
-    _source_terms.push_back(std::move(per_process_source_terms));
+    auto& per_process_sts = _source_term_collections[process_id];
+    per_process_sts.addSourceTermsForProcessVariables(
+        per_process_variables, dof_table, _integration_order);
 }
 
 void Process::initializeBoundaryConditions()
@@ -188,11 +192,7 @@ void Process::assemble(const double t, GlobalVector const& x, GlobalMatrix& M,
         (_coupled_solutions) != nullptr ? _coupled_solutions->process_id : 0;
     _boundary_conditions[pcs_id].applyNaturalBC(t, x, K, b);
 
-    auto& source_terms_per_pcs = _source_terms[pcs_id];
-    for (auto& st : source_terms_per_pcs)
-    {
-        st->integrateNodalSourceTerm(t, b);
-    }
+    _source_term_collections[pcs_id].integrateNodalSourceTerms(t, b);
 }
 
 void Process::assembleWithJacobian(const double t, GlobalVector const& x,

diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
@@ -16,6 +16,7 @@
 #include "NumLib/ODESolver/ODESystem.h"
 #include "NumLib/ODESolver/TimeDiscretization.h"
 #include "ProcessLib/BoundaryCondition/BoundaryConditionCollection.h"
+#include "ProcessLib/SourceTerms/SourceTermCollection.h"
 #include "ProcessLib/Output/CachedSecondaryVariable.h"
 #include "ProcessLib/Output/ExtrapolatorData.h"
 #include "ProcessLib/Output/SecondaryVariable.h"
@@ -278,9 +279,10 @@ class Process
     /// scheme, the size of vector is the number of the coupled processes.
     std::vector<BoundaryConditionCollection> _boundary_conditions;
 
-    /// Vector for nodal source terms. The outer vector is for processes,
-    /// which has the same size as that for boundary conditions.
-    std::vector<std::vector<std::unique_ptr<NodalSourceTerm>>> _source_terms;
+    /// Vector for nodal source term collections. For the monolithic scheme
+    /// or a single process, the size of the vector is one. For the staggered
+    /// scheme, the size of vector is the number of the coupled processes.
+    std::vector<SourceTermCollection> _source_term_collections;
 
     ExtrapolatorData _extrapolator_data;
 };

diff --git a/ProcessLib/ProcessVariable.cpp b/ProcessLib/ProcessVariable.cpp
@@ -178,14 +178,15 @@ std::vector<std::unique_ptr<NodalSourceTerm>>
 ProcessVariable::createSourceTerms(
     const NumLib::LocalToGlobalIndexMap& dof_table,
     const int variable_id,
-    unsigned const integration_order)
+    unsigned const integration_order,
+    std::vector<std::unique_ptr<ParameterBase>> const& parameters)
 {
     std::vector<std::unique_ptr<NodalSourceTerm>> source_terms;
 
     for (auto& config : _source_term_configs)
         source_terms.emplace_back(_source_term_builder->createSourceTerm(
             config, dof_table, _mesh, variable_id, integration_order,
-            _shapefunction_order));
+            _shapefunction_order, parameters));
 
     return source_terms;
 }

diff --git a/ProcessLib/ProcessVariable.h b/ProcessLib/ProcessVariable.h
@@ -57,7 +57,8 @@ class ProcessVariable
 
     std::vector<std::unique_ptr<NodalSourceTerm>> createSourceTerms(
         const NumLib::LocalToGlobalIndexMap& dof_table, const int variable_id,
-        unsigned const integration_order);
+        unsigned const integration_order,
+        std::vector<std::unique_ptr<ParameterBase>> const& parameters);
 
     Parameter<double> const& getInitialCondition() const
     {

diff --git a/ProcessLib/SourceTerms/CreateNodalSourceTerm.cpp b/ProcessLib/SourceTerms/CreateNodalSourceTerm.cpp
@@ -18,18 +18,21 @@ namespace ProcessLib
 std::unique_ptr<NodalSourceTerm> createNodalSourceTerm(
     BaseLib::ConfigTree const& config,
     const NumLib::LocalToGlobalIndexMap& dof_table, std::size_t const mesh_id,
-    std::size_t const node_id, const int variable_id, const int component_id)
+    std::size_t const node_id, const int variable_id, const int component_id,
+    std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const& parameters)
 {
     DBUG("Constructing NodalSourceTerm from config.");
     //! \ogs_file_param{prj__process_variables__process_variable__source_terms__source_term__type}
     config.checkConfigParameter("type", "Nodal");
 
-    //! \ogs_file_param{prj__process_variables__process_variable__source_terms__source_term__Nodal__value}
-    auto const nodal_value = config.getConfigParameter<double>("value");
-    DBUG("Using value %f as nodal source term", nodal_value);
+    //! \ogs_file_param{prj__process_variables__process_variable__source_terms__source_term__Nodal__parameter}
+    auto const param_name = config.getConfigParameter<std::string>("parameter");
+    DBUG("Using parameter %s as nodal source term.", param_name.c_str());
+
+    auto& param = findParameter<double>(param_name, parameters, 1);
 
     return std::make_unique<NodalSourceTerm>(
-        dof_table, mesh_id, node_id, variable_id, component_id, nodal_value);
+        dof_table, mesh_id, node_id, variable_id, component_id, param);
 }
 
 }  // namespace ProcessLib
diff --git a/ProcessLib/SourceTerms/CreateNodalSourceTerm.h b/ProcessLib/SourceTerms/CreateNodalSourceTerm.h
@@ -17,6 +17,7 @@ namespace ProcessLib
 std::unique_ptr<NodalSourceTerm> createNodalSourceTerm(
     BaseLib::ConfigTree const& config,
     const NumLib::LocalToGlobalIndexMap& dof_table, std::size_t mesh_id,
-    std::size_t const node_id, const int variable_id, const int component_id);
+    std::size_t const node_id, const int variable_id, const int component_id,
+    std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const& parameters);
 
 }   // namespace ProcessLib
diff --git a/ProcessLib/SourceTerms/NodalSourceTerm.cpp b/ProcessLib/SourceTerms/NodalSourceTerm.cpp
@@ -17,26 +17,30 @@ NodalSourceTerm::NodalSourceTerm(const NumLib::LocalToGlobalIndexMap& dof_table,
                                  std::size_t const mesh_id,
                                  std::size_t const node_id,
                                  const int variable_id, const int component_id,
-                                 double value)
+                                 Parameter<double> const& parameter)
     : _dof_table(dof_table),
       _mesh_id(mesh_id),
       _node_id(node_id),
       _variable_id(variable_id),
       _component_id(component_id),
-      _value(value)
+      _parameter(parameter)
 {
     DBUG("Create NodalSourceTerm.");
 }
 
-void NodalSourceTerm::integrateNodalSourceTerm(const double /*t*/,
+void NodalSourceTerm::integrateNodalSourceTerm(const double t,
                                                GlobalVector& b) const
 {
     DBUG("Assemble NodalSourceTerm.");
 
     MeshLib::Location const l{_mesh_id, MeshLib::MeshItemType::Node, _node_id};
     auto const index =
         _dof_table.getGlobalIndex(l, _variable_id, _component_id);
-    b.add(index, _value);
+
+    SpatialPosition pos;
+    pos.setNodeID(_node_id);
+
+    b.add(index, _parameter(t, pos).front());
 }
 
 }  // namespace ProcessLib
diff --git a/ProcessLib/SourceTerms/NodalSourceTerm.h b/ProcessLib/SourceTerms/NodalSourceTerm.h
@@ -10,6 +10,7 @@
 #pragma once
 
 #include "NumLib/DOF/LocalToGlobalIndexMap.h"
+#include "ProcessLib/Parameter/Parameter.h"
 
 namespace ProcessLib
 {
@@ -19,7 +20,7 @@ class NodalSourceTerm final
     NodalSourceTerm(const NumLib::LocalToGlobalIndexMap& dof_table,
                     std::size_t const mesh_id, std::size_t const node_id,
                     const int variable_id, const int component_id,
-                    double value);
+                    Parameter<double> const& parameter);
 
     void integrateNodalSourceTerm(
         const double t,
@@ -31,7 +32,7 @@ class NodalSourceTerm final
     std::size_t const _node_id;
     int const _variable_id;
     int const _component_id;
-    double const _value;
+    Parameter<double> const& _parameter;
 };
 
 }   // namespace ProcessLib
diff --git a/ProcessLib/SourceTerms/SourceTermBuilder.cpp b/ProcessLib/SourceTerms/SourceTermBuilder.cpp
@@ -22,15 +22,17 @@ std::unique_ptr<NodalSourceTerm> SourceTermBuilder::createSourceTerm(
     const SourceTermConfig& config,
     const NumLib::LocalToGlobalIndexMap& dof_table, const MeshLib::Mesh& mesh,
     const int variable_id, const unsigned integration_order,
-    const unsigned shapefunction_order)
+    const unsigned shapefunction_order,
+    std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const& parameters)
 {
     //! \ogs_file_param{prj__process_variables__process_variable__source_terms__source_term__type}
     auto const type = config.config.peekConfigParameter<std::string>("type");
 
     if (type == "Nodal")
     {
         return createNodalSourceTerm(config, dof_table, mesh, variable_id,
-                                     integration_order, shapefunction_order);
+                                     integration_order, shapefunction_order,
+                                     parameters);
     }
 
     OGS_FATAL("Unknown source term type: `%s'.", type.c_str());
@@ -40,7 +42,8 @@ std::unique_ptr<NodalSourceTerm> SourceTermBuilder::createNodalSourceTerm(
     const SourceTermConfig& config,
     const NumLib::LocalToGlobalIndexMap& dof_table, const MeshLib::Mesh& mesh,
     const int variable_id, const unsigned /*integration_order*/,
-    const unsigned /*shapefunction_order*/)
+    const unsigned /*shapefunction_order*/,
+    std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const& parameters)
 {
     std::unique_ptr<MeshGeoToolsLib::SearchLength> search_length_algorithm =
         MeshGeoToolsLib::createSearchLengthAlgorithm(config.config, mesh);
@@ -88,11 +91,11 @@ std::unique_ptr<NodalSourceTerm> SourceTermBuilder::createNodalSourceTerm(
     if (ids.size() != 1)
         OGS_FATAL(
             "Found %d nodes for nodal source term, but exactly one node is "
-            "required.");
+            "required.", ids.size());
 
     return ProcessLib::createNodalSourceTerm(
         config.config, dof_table, mesh.getID(), ids[0], variable_id,
-        *config.component_id);
+        *config.component_id, parameters);
 }
 
 }  // ProcessLib
diff --git a/ProcessLib/SourceTerms/SourceTermBuilder.h b/ProcessLib/SourceTerms/SourceTermBuilder.h
@@ -10,6 +10,7 @@
 #pragma once
 
 #include "NodalSourceTerm.h"
+#include "ProcessLib/Parameter/Parameter.h"
 
 namespace GeoLib
 {
@@ -47,14 +48,18 @@ class SourceTermBuilder
         const SourceTermConfig& config,
         const NumLib::LocalToGlobalIndexMap& dof_table,
         const MeshLib::Mesh& mesh, const int variable_id,
-        const unsigned integration_order, const unsigned shapefunction_order);
+        const unsigned integration_order, const unsigned shapefunction_order,
+        std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const&
+            parameters);
 
 protected:
     virtual std::unique_ptr<NodalSourceTerm> createNodalSourceTerm(
         const SourceTermConfig& config,
         const NumLib::LocalToGlobalIndexMap& dof_table,
         const MeshLib::Mesh& mesh, const int variable_id,
-        const unsigned integration_order, const unsigned shapefunction_order);
+        const unsigned integration_order, const unsigned shapefunction_order,
+        std::vector<std::unique_ptr<ProcessLib::ParameterBase>> const&
+            parameters);
 };
 
 }  // ProcessLib
diff --git a/ProcessLib/SourceTerms/SourceTermCollection.cpp b/ProcessLib/SourceTerms/SourceTermCollection.cpp
@@ -0,0 +1,39 @@
+/**
+ * \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 "SourceTermCollection.h"
+
+namespace ProcessLib
+{
+void SourceTermCollection::addSourceTermsForProcessVariables(
+    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 sts = pv.createSourceTerms(dof_table, variable_id,
+                                        integration_order, _parameters);
+
+        std::move(sts.begin(), sts.end(), std::back_inserter(_source_terms));
+    }
+}
+
+void SourceTermCollection::integrateNodalSourceTerms(const double t,
+                                                     GlobalVector& b) const
+{
+    for (auto const& st : _source_terms)
+        st->integrateNodalSourceTerm(t, b);
+}
+
+}
diff --git a/ProcessLib/SourceTerms/SourceTermCollection.h b/ProcessLib/SourceTerms/SourceTermCollection.h
@@ -0,0 +1,39 @@
+/**
+ * \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 "ProcessLib/ProcessVariable.h"
+#include "ProcessLib/SourceTerms/NodalSourceTerm.h"
+
+namespace ProcessLib
+{
+class SourceTermCollection final
+{
+public:
+    SourceTermCollection(
+        std::vector<std::unique_ptr<ParameterBase>> const& parameters)
+        : _parameters(parameters)
+    {
+    }
+
+    void integrateNodalSourceTerms(const double t, GlobalVector& b) const;
+
+    void addSourceTermsForProcessVariables(
+        std::vector<std::reference_wrapper<ProcessVariable>> const&
+            process_variables,
+        NumLib::LocalToGlobalIndexMap const& dof_table,
+        unsigned const integration_order);
+
+private:
+    std::vector<std::unique_ptr<NodalSourceTerm>> _source_terms;
+    std::vector<std::unique_ptr<ParameterBase>> const& _parameters;
+};
+
+}  // ProcessLib
diff --git a/Tests/Data/Elliptic/circle_radius_1/circle_1e1_axi.prj b/Tests/Data/Elliptic/circle_radius_1/circle_1e1_axi.prj
@@ -84,6 +84,11 @@
             <type>Constant</type>
             <value>0</value>
         </parameter>
+        <parameter>
+            <name>pressure_source_term</name>
+            <type>Constant</type>
+            <value>1</value>
+        </parameter>
     </parameters>
     <process_variables>
         <process_variable>
@@ -104,7 +109,7 @@
                     <geometrical_set>geometry</geometrical_set>
                     <geometry>inner</geometry>
                     <type>Nodal</type>
-                    <value>1</value>
+                    <parameter>pressure_source_term</parameter>
                 </source_term>
             </source_terms>
         </process_variable>

diff --git a/Tests/Data/Elliptic/circle_radius_1/circle_1e2_axi.prj b/Tests/Data/Elliptic/circle_radius_1/circle_1e2_axi.prj
@@ -84,6 +84,11 @@
             <type>Constant</type>
             <value>0</value>
         </parameter>
+        <parameter>
+            <name>pressure_source_term</name>
+            <type>Constant</type>
+            <value>1</value>
+        </parameter>
     </parameters>
     <process_variables>
         <process_variable>
@@ -104,7 +109,7 @@
                     <geometrical_set>geometry</geometrical_set>
                     <geometry>inner</geometry>
                     <type>Nodal</type>
-                    <value>1</value>
+                    <parameter>pressure_source_term</parameter>
                 </source_term>
             </source_terms>
         </process_variable>