diff --git a/Applications/ApplicationsLib/ProjectData.cpp b/Applications/ApplicationsLib/ProjectData.cpp
index d6bb86d88cc0dc1aab94585cebb93e66dda52c68..a035a3fada6418ef29c00aaed21f36225165553d 100644
--- a/Applications/ApplicationsLib/ProjectData.cpp
+++ b/Applications/ApplicationsLib/ProjectData.cpp
@@ -36,6 +36,7 @@
#include "ProcessLib/HeatConduction/CreateHeatConductionProcess.h"
#include "ProcessLib/RichardsFlow/CreateRichardsFlowProcess.h"
#include "ProcessLib/SmallDeformation/CreateSmallDeformationProcess.h"
+#include "ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.h"
#include "ProcessLib/TES/CreateTESProcess.h"
namespace detail
@@ -330,6 +331,23 @@ void ProjectData::parseProcesses(BaseLib::ConfigTree const& processes_config,
"given dimension");
}
}
+ else if (type == "SMALL_DEFORMATION_WITH_LIE")
+ {
+ switch (process_config.getConfigParameter<int>("dimension"))
+ {
+ case 2:
+ process = ProcessLib::SmallDeformationWithLIE::
+ createSmallDeformationProcess<2>(
+ *_mesh_vec[0], std::move(jacobian_assembler),
+ _process_variables, _parameters, integration_order,
+ process_config);
+ break;
+ default:
+ OGS_FATAL(
+ "SMALL_DEFORMATION_WITH_LIE process does not support "
+ "given dimension");
+ }
+ }
else if (type == "RICHARDS_FLOW")
{
process = ProcessLib::RichardsFlow::createRichardsFlowProcess(
diff --git a/Applications/CLI/Tests.cmake b/Applications/CLI/Tests.cmake
index f50ef2788a000d97657b4769ebe0ab8f4ebe5215..424b9b6f0c43d7b9328c41e6c040ac72174739a8 100644
--- a/Applications/CLI/Tests.cmake
+++ b/Applications/CLI/Tests.cmake
@@ -503,6 +503,20 @@ if(NOT OGS_USE_MPI)
# EXECUTABLE_ARGS tes-inert-wedge.prj
# )
+ # LIE; Small deformation
+ AddTest(
+ NAME LIE_M_single_joint
+ PATH LIE/Mechanics
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS single_joint.prj
+ WRAPPER time
+ TESTER vtkdiff
+ ABSTOL 1e-16 RELTOL 1e-16
+ DIFF_DATA
+ single_joint_expected_pcs_0_ts_1_t_1.000000.vtu single_joint_pcs_0_ts_1_t_1.000000.vtu displacement displacement
+ single_joint_expected_pcs_0_ts_1_t_1.000000.vtu single_joint_pcs_0_ts_1_t_1.000000.vtu displacement_jump1 displacement_jump1
+ )
+
else()
# MPI groundwater flow tests
AddTest(
diff --git a/NumLib/DOF/LocalToGlobalIndexMap.cpp b/NumLib/DOF/LocalToGlobalIndexMap.cpp
index f7dd2c1625ae547ea33433b83a66de39000d9c3e..7c046a9a35d49ed423d5204cf55d01b5b0f4f134 100644
--- a/NumLib/DOF/LocalToGlobalIndexMap.cpp
+++ b/NumLib/DOF/LocalToGlobalIndexMap.cpp
@@ -313,6 +313,7 @@ std::ostream& operator<<(std::ostream& os, LocalToGlobalIndexMap const& map)
<< " rows\n";
for (std::size_t e=0; e<map.size(); ++e)
{
+ os << "== e " << e << " ==\n";
for (std::size_t c=0; c<map.getNumberOfComponents(); ++c)
{
auto const& line = map._rows(e, c);
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index 4615d204d3697d2f80b746b3280ce50b56351bdb..7ee405f8d4ac0456d12ab1b3e3dd1e42e47bc539 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -154,8 +154,10 @@ private:
return NumLib::IterationResult::SUCCESS;
}
- void constructDofTable();
+protected:
+ virtual void constructDofTable();
+private:
void initializeExtrapolator();
/// Finishes the \c _named_function_caller and adds a secondary variable for
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess-fwd.h b/ProcessLib/SmallDeformation/SmallDeformationProcess-fwd.h
index 0cb7eaabdc8f4b7777d03d110b0f27210f00d4c1..d5110cae60699164e675e5e62eb2d202303fcb77 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess-fwd.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess-fwd.h
@@ -1,6 +1,6 @@
/**
* \copyright
- * Copyright (c) 2012-2015, OpenGeoSys Community (http://www.opengeosys.org)
+ * 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
diff --git a/ProcessLib/SmallDeformationWithLIE/Common/FractureProperty.h b/ProcessLib/SmallDeformationWithLIE/Common/FractureProperty.h
index e6295eed80e1f297da137c0ed08989ce305c6066..3f8fff59c60ce341fcd93d3320db4533f4cc70e3 100644
--- a/ProcessLib/SmallDeformationWithLIE/Common/FractureProperty.h
+++ b/ProcessLib/SmallDeformationWithLIE/Common/FractureProperty.h
@@ -17,7 +17,6 @@
#include "ProcessLib/Parameter/Parameter.h"
-#include "FractureProperty.h"
#include "Utils.h"
diff --git a/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.cpp b/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..1b56136cc1dca4163fe89daa685791e3910c03ea
--- /dev/null
+++ b/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.cpp
@@ -0,0 +1,177 @@
+/**
+ * \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 "CreateSmallDeformationProcess.h"
+
+#include <cassert>
+
+#include "MaterialLib/SolidModels/CreateLinearElasticIsotropic.h"
+#include "MaterialLib/FractureModels/CreateLinearElasticIsotropic.h"
+#include "MaterialLib/FractureModels/CreateMohrCoulomb.h"
+
+#include "ProcessLib/Utils/ParseSecondaryVariables.h"
+#include "ProcessLib/Utils/ProcessUtils.h" // required for findParameter
+
+#include "SmallDeformationProcess.h"
+#include "SmallDeformationProcessData.h"
+
+namespace ProcessLib
+{
+namespace SmallDeformationWithLIE
+{
+template <int DisplacementDim>
+class SmallDeformationProcess;
+
+extern template class SmallDeformationProcess<2>;
+
+template <int DisplacementDim>
+std::unique_ptr<Process>
+createSmallDeformationProcess(
+ MeshLib::Mesh& mesh,
+ std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
+ std::vector<ProcessVariable> const& variables,
+ std::vector<std::unique_ptr<ParameterBase>> const& parameters,
+ unsigned const integration_order,
+ BaseLib::ConfigTree const& config)
+{
+ //! \ogs_file_param{process__type}
+ config.checkConfigParameter("type", "SMALL_DEFORMATION_WITH_LIE");
+ DBUG("Create SmallDeformationProcess with LIE.");
+
+ // Process variables
+ auto const pv_conf = config.getConfigSubtree("process_variables");
+ auto range = pv_conf.getConfigParameterList<std::string>("process_variable");
+ std::vector<std::reference_wrapper<ProcessVariable>> process_variables;
+ for (std::string const& pv_name : range)
+ {
+ if (pv_name != "displacement"
+ && pv_name.find("displacement_jump")==std::string::npos)
+ OGS_FATAL("Found a process variable name '%s'. It should be 'displacement' or 'displacement_jumpN'");
+ auto variable = std::find_if(
+ variables.cbegin(), variables.cend(),
+ [&pv_name](ProcessVariable const& v) { return v.getName() == pv_name; });
+
+ if (variable == variables.end())
+ {
+ OGS_FATAL(
+ "Could not find process variable '%s' in the provided variables "
+ "list for config tag <%s>.",
+ pv_name.c_str(), "process_variable");
+ }
+ DBUG("Found process variable \'%s\' for config tag <%s>.",
+ variable->getName().c_str(), "process_variable");
+
+ process_variables.emplace_back(const_cast<ProcessVariable&>(*variable));
+ }
+ if (process_variables.size() > 2)
+ OGS_FATAL("Currently only one displacement jump is supported");
+
+ DBUG("Associate displacement with process variable \'%s\'.",
+ process_variables.back().get().getName().c_str());
+
+ if (process_variables.back().get().getNumberOfComponents() !=
+ DisplacementDim)
+ {
+ OGS_FATAL(
+ "Number of components of the process variable '%s' is different "
+ "from the displacement dimension: got %d, expected %d",
+ process_variables.back().get().getName().c_str(),
+ process_variables.back().get().getNumberOfComponents(),
+ DisplacementDim);
+ }
+
+ // Constitutive relation.
+ // read type;
+ auto const constitutive_relation_config =
+ //! \ogs_file_param{process__SMALL_DEFORMATION_WITH_LIE__constitutive_relation}
+ config.getConfigSubtree("constitutive_relation");
+
+ auto const type =
+ constitutive_relation_config.peekConfigParameter<std::string>("type");
+
+ std::unique_ptr<MaterialLib::Solids::MechanicsBase<DisplacementDim>> material = nullptr;
+ if (type == "LinearElasticIsotropic")
+ {
+ material = MaterialLib::Solids::createLinearElasticIsotropic<DisplacementDim>(
+ parameters, constitutive_relation_config);
+ }
+ else
+ {
+ OGS_FATAL(
+ "Cannot construct constitutive relation of given type \'%s\'.",
+ type.c_str());
+ }
+
+ // Fracture constitutive relation.
+ // read type;
+ auto const fracture_constitutive_relation_config =
+ //! \ogs_file_param{process__SMALL_DEFORMATION_WITH_LIE__constitutive_relation}
+ config.getConfigSubtree("fracture_constitutive_relation");
+
+ auto const frac_type =
+ fracture_constitutive_relation_config.peekConfigParameter<std::string>("type");
+
+ std::unique_ptr<MaterialLib::Fracture::FractureModelBase<DisplacementDim>> fracture_model = nullptr;
+ if (frac_type == "LinearElasticIsotropic")
+ {
+ fracture_model = MaterialLib::Fracture::createLinearElasticIsotropic<DisplacementDim>(
+ parameters, fracture_constitutive_relation_config);
+ }
+ else if (frac_type == "MohrCoulomb")
+ {
+ fracture_model = MaterialLib::Fracture::createMohrCoulomb<DisplacementDim>(
+ parameters, fracture_constitutive_relation_config);
+ }
+ else
+ {
+ OGS_FATAL(
+ "Cannot construct fracture constitutive relation of given type \'%s\'.",
+ frac_type.c_str());
+ }
+
+ // Fracture properties
+ //! \ogs_file_param{process__SMALL_DEFORMATION_WITH_LIE__fracture_properties}
+ auto fracture_properties_config = config.getConfigSubtree("fracture_properties");
+ auto ¶_b0 = ProcessLib::findParameter<double>(fracture_properties_config, "initial_aperture", parameters, 1);
+ std::unique_ptr<FractureProperty> frac_prop(new FractureProperty());
+ frac_prop->mat_id = fracture_properties_config.getConfigParameter<int>("material_id");
+ frac_prop->aperture0 = ¶_b0;
+
+
+ SmallDeformationProcessData<DisplacementDim> process_data(
+ std::move(material), std::move(fracture_model), std::move(frac_prop));
+
+ SecondaryVariableCollection secondary_variables;
+
+ NumLib::NamedFunctionCaller named_function_caller(
+ {"SmallDeformation_displacement"});
+
+ ProcessLib::parseSecondaryVariables(config, secondary_variables,
+ named_function_caller);
+
+ return std::unique_ptr<SmallDeformationProcess<DisplacementDim>>{
+ new SmallDeformationProcess<DisplacementDim>{
+ mesh, std::move(jacobian_assembler), parameters, integration_order,
+ std::move(process_variables), std::move(process_data),
+ std::move(secondary_variables), std::move(named_function_caller)}};
+}
+
+template
+std::unique_ptr<Process>
+createSmallDeformationProcess<2>(
+ MeshLib::Mesh& mesh,
+ std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
+ std::vector<ProcessVariable> const& variables,
+ std::vector<std::unique_ptr<ParameterBase>> const& parameters,
+ unsigned const integration_order,
+ BaseLib::ConfigTree const& config);
+
+} // namespace SmallDeformationWithLIE
+} // namespace ProcessLib
+
diff --git a/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.h b/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.h
new file mode 100644
index 0000000000000000000000000000000000000000..0b30c314b38264fa1fb6d15628c39ea2dd10d3e9
--- /dev/null
+++ b/ProcessLib/SmallDeformationWithLIE/CreateSmallDeformationProcess.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_SMALLDEFORMATION_WITH_LIE_CREATESMALLDEFORMATIONPROCESS_H_
+#define PROCESS_LIB_SMALLDEFORMATION_WITH_LIE_CREATESMALLDEFORMATIONPROCESS_H_
+
+#include <memory>
+#include <vector>
+
+#include "ProcessLib/Process.h"
+
+namespace ProcessLib
+{
+namespace SmallDeformationWithLIE
+{
+template <int DisplacementDim>
+std::unique_ptr<Process> createSmallDeformationProcess(
+ MeshLib::Mesh& mesh,
+ std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
+ std::vector<ProcessVariable> const& variables,
+ std::vector<std::unique_ptr<ParameterBase>> const& parameters,
+ unsigned const integration_order,
+ BaseLib::ConfigTree const& config);
+
+} // namespace SmallDeformationWithLIE
+} // namespace ProcessLib
+
+#endif // PROCESS_LIB_SMALLDEFORMATION_WITH_LIE_CREATESMALLDEFORMATIONPROCESS_H_
diff --git a/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h b/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
index a1f6f349fa2fdd4e6cddbf4b09e9f4c6a2f28162..3a5a0dd878f3b94c195f53091f24b6678945c2e4 100644
--- a/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
+++ b/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
@@ -147,7 +147,10 @@ private:
cache.reserve(_ip_data.size());
for (auto const& ip_data : _ip_data) {
- cache.push_back(ip_data._sigma[component]);
+ if (component < 3) // xx, yy, zz components
+ cache.push_back(ip_data._sigma[component]);
+ else // mixed xy, yz, xz components
+ cache.push_back(ip_data._sigma[component] / std::sqrt(2));
}
return cache;
diff --git a/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h b/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
index 3c22413f85556c47f6dcbdad13fef189c531ecd5..d3356ea8693336e78bbd5dde4cd3b5eb964a6ea9 100644
--- a/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
+++ b/ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
@@ -148,7 +148,10 @@ private:
cache.reserve(_ip_data.size());
for (auto const& ip_data : _ip_data) {
- cache.push_back(ip_data._sigma[component]);
+ if (component < 3) // xx, yy, zz components
+ cache.push_back(ip_data._sigma[component]);
+ else // mixed xy, yz, xz components
+ cache.push_back(ip_data._sigma[component] / std::sqrt(2));
}
return cache;
diff --git a/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess-fwd.h b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess-fwd.h
new file mode 100644
index 0000000000000000000000000000000000000000..58ac118ebcca462ff21959e35d8239b3580b25f2
--- /dev/null
+++ b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess-fwd.h
@@ -0,0 +1,17 @@
+/**
+ * \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_SMALLDEFORMATION_WITH_LIE_PROCESS_FWD_H_
+#define PROCESSLIB_SMALLDEFORMATION_WITH_LIE_PROCESS_FWD_H_
+
+#include "SmallDeformationProcess.h"
+
+extern template class ProcessLib::SmallDeformationWithLIE::SmallDeformationProcess<2>;
+
+#endif // PROCESSLIB_SMALLDEFORMATION_WITH_LIE_PROCESS_FWD_H_
diff --git a/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..5efec3fc6cc60b47c4cd2b71b526b32723d44274
--- /dev/null
+++ b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.cpp
@@ -0,0 +1,278 @@
+/**
+ * \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 "SmallDeformationProcess-fwd.h"
+#include "SmallDeformationProcess.h"
+
+#include <algorithm>
+#include <cassert>
+#include <iostream>
+#include <vector>
+
+#include "MeshLib/Elements/Element.h"
+#include "MeshLib/ElementCoordinatesMappingLocal.h"
+#include "MeshLib/Mesh.h"
+#include "MeshLib/Node.h"
+#include "MeshLib/Properties.h"
+
+#include "NumLib/DOF/LocalToGlobalIndexMap.h"
+
+#include "ProcessLib/SmallDeformationWithLIE/BoundaryCondition/BoundaryConditionBuilder.h"
+#include "ProcessLib/SmallDeformationWithLIE/Common/LevelSetFunction.h"
+#include "ProcessLib/SmallDeformationWithLIE/Common/MeshUtils.h"
+#include "ProcessLib/SmallDeformationWithLIE/Common/Utils.h"
+#include "ProcessLib/SmallDeformationWithLIE/LocalAssembler/LocalAssemblerDataMatrix.h"
+#include "ProcessLib/SmallDeformationWithLIE/LocalAssembler/LocalAssemblerDataMatrixNearFracture.h"
+#include "ProcessLib/SmallDeformationWithLIE/LocalAssembler/LocalAssemblerDataFracture.h"
+
+namespace ProcessLib
+{
+namespace SmallDeformationWithLIE
+{
+
+template <int DisplacementDim>
+SmallDeformationProcess<DisplacementDim>::SmallDeformationProcess(
+ MeshLib::Mesh& mesh,
+ std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&&
+ jacobian_assembler,
+ std::vector<std::unique_ptr<ParameterBase>> const& parameters,
+ unsigned const integration_order,
+ std::vector<std::reference_wrapper<ProcessVariable>>&&
+ process_variables,
+ SmallDeformationProcessData<DisplacementDim>&& process_data,
+ SecondaryVariableCollection&& secondary_variables,
+ NumLib::NamedFunctionCaller&& named_function_caller)
+ : Process(mesh, std::move(jacobian_assembler), parameters,
+ integration_order,
+ std::move(process_variables),
+ std::move(secondary_variables),
+ std::move(named_function_caller)),
+ _process_data(std::move(process_data))
+{
+ getFractureMatrixDataInMesh(mesh,
+ _vec_matrix_elements,
+ _vec_fracture_elements,
+ _vec_fracture_matrix_elements,
+ _vec_fracture_nodes);
+
+ // set fracture property assuming a fracture forms a straight line
+ setFractureProperty(DisplacementDim,
+ *_vec_fracture_elements[0],
+ *_process_data._fracture_property.get());
+
+ // need to use a custom Neumann BC assembler for displacement jumps
+ for (ProcessVariable& pv : getProcessVariables())
+ {
+ if (pv.getName().find("displacement_jump") == std::string::npos)
+ continue;
+ pv.setBoundaryConditionBuilder(
+ std::unique_ptr<ProcessLib::BoundaryConditionBuilder>(
+ new BoundaryConditionBuilder(*_process_data._fracture_property.get())));
+ }
+
+}
+
+
+template <int DisplacementDim>
+void SmallDeformationProcess<DisplacementDim>::constructDofTable()
+{
+ //------------------------------------------------------------
+ // prepare mesh subsets to define DoFs
+ //------------------------------------------------------------
+ // for extrapolation
+ _mesh_subset_all_nodes.reset(new MeshLib::MeshSubset(_mesh, &_mesh.getNodes()));
+ // regular u
+ _mesh_subset_matrix_nodes.reset(new MeshLib::MeshSubset(_mesh, &_mesh.getNodes()));
+ // u jump
+ _mesh_subset_fracture_nodes.reset(new MeshLib::MeshSubset(_mesh, &_vec_fracture_nodes));
+
+ // Collect the mesh subsets in a vector.
+ std::vector<std::unique_ptr<MeshLib::MeshSubsets>> all_mesh_subsets;
+ std::generate_n(
+ std::back_inserter(all_mesh_subsets),
+ DisplacementDim,
+ [&]() {
+ return std::unique_ptr<MeshLib::MeshSubsets>{
+ new MeshLib::MeshSubsets{_mesh_subset_matrix_nodes.get()}};
+ });
+ std::generate_n(
+ std::back_inserter(all_mesh_subsets),
+ DisplacementDim,
+ [&]() {
+ return std::unique_ptr<MeshLib::MeshSubsets>{
+ new MeshLib::MeshSubsets{_mesh_subset_fracture_nodes.get()}};
+ });
+
+ std::vector<unsigned> const vec_n_components(2, DisplacementDim);
+
+ std::vector<std::vector<MeshLib::Element*>const*> vec_var_elements;
+ vec_var_elements.push_back(&_vec_matrix_elements);
+ vec_var_elements.push_back(&_vec_fracture_matrix_elements);
+
+ _local_to_global_index_map.reset(
+ new NumLib::LocalToGlobalIndexMap(
+ std::move(all_mesh_subsets),
+ vec_n_components,
+ vec_var_elements,
+ NumLib::ComponentOrder::BY_COMPONENT));
+}
+
+
+template <int DisplacementDim>
+void SmallDeformationProcess<DisplacementDim>::initializeConcreteProcess(
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ MeshLib::Mesh const& mesh,
+ unsigned const integration_order)
+{
+ ProcessLib::SmallDeformationWithLIE::createLocalAssemblers
+ <DisplacementDim,
+ LocalAssemblerDataMatrix,
+ LocalAssemblerDataMatrixNearFracture,
+ LocalAssemblerDataFracture>(
+ mesh.getDimension(), mesh.getElements(), dof_table,
+ _local_assemblers, mesh.isAxiallySymmetric(), integration_order,
+ _process_data);
+
+ // TODO move the two data members somewhere else.
+ // for extrapolation of secondary variables
+ std::vector<std::unique_ptr<MeshLib::MeshSubsets>>
+ all_mesh_subsets_single_component;
+ all_mesh_subsets_single_component.emplace_back(
+ new MeshLib::MeshSubsets(_mesh_subset_all_nodes.get()));
+ _local_to_global_index_map_single_component.reset(
+ new NumLib::LocalToGlobalIndexMap(
+ std::move(all_mesh_subsets_single_component),
+ // by location order is needed for output
+ NumLib::ComponentOrder::BY_LOCATION));
+
+ Base::_secondary_variables.addSecondaryVariable(
+ "sigma_xx", 1,
+ makeExtrapolator(
+ getExtrapolator(), _local_assemblers,
+ &SmallDeformationLocalAssemblerInterface::getIntPtSigmaXX));
+
+ Base::_secondary_variables.addSecondaryVariable(
+ "sigma_yy", 1,
+ makeExtrapolator(
+ getExtrapolator(), _local_assemblers,
+ &SmallDeformationLocalAssemblerInterface::getIntPtSigmaYY));
+
+ Base::_secondary_variables.addSecondaryVariable(
+ "sigma_zz", 1,
+ makeExtrapolator(
+ getExtrapolator(), _local_assemblers,
+ &SmallDeformationLocalAssemblerInterface::getIntPtSigmaZZ));
+
+ Base::_secondary_variables.addSecondaryVariable(
+ "sigma_xy", 1,
+ makeExtrapolator(
+ getExtrapolator(), _local_assemblers,
+ &SmallDeformationLocalAssemblerInterface::getIntPtSigmaXY));
+
+ auto mesh_prop_sigma_xx = const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "stress_xx", MeshLib::MeshItemType::Cell);
+ mesh_prop_sigma_xx->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_stress_xx = mesh_prop_sigma_xx;
+
+ auto mesh_prop_sigma_yy = const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "stress_yy", MeshLib::MeshItemType::Cell);
+ mesh_prop_sigma_yy->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_stress_yy = mesh_prop_sigma_yy;
+
+ auto mesh_prop_sigma_xy = const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "stress_xy", MeshLib::MeshItemType::Cell);
+ mesh_prop_sigma_xy->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_stress_xy = mesh_prop_sigma_xy;
+
+ auto mesh_prop_epsilon_xx =
+ const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "strain_xx", MeshLib::MeshItemType::Cell);
+ mesh_prop_epsilon_xx->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_strain_xx = mesh_prop_epsilon_xx;
+
+ auto mesh_prop_epsilon_yy =
+ const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "strain_yy", MeshLib::MeshItemType::Cell);
+ mesh_prop_epsilon_yy->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_strain_yy = mesh_prop_epsilon_yy;
+
+ auto mesh_prop_epsilon_xy =
+ const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "strain_xy", MeshLib::MeshItemType::Cell);
+ mesh_prop_epsilon_xy->resize(mesh.getNumberOfElements());
+ _process_data._mesh_prop_strain_xy = mesh_prop_epsilon_xy;
+
+ auto mesh_prop_levelset = const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "levelset1", MeshLib::MeshItemType::Cell);
+ mesh_prop_levelset->resize(mesh.getNumberOfElements());
+ for (MeshLib::Element const* e : _mesh.getElements())
+ {
+ if (e->getDimension() < DisplacementDim)
+ continue;
+
+ double const levelsets =
+ calculateLevelSetFunction(*_process_data._fracture_property,
+ e->getCenterOfGravity().getCoords());
+ (*mesh_prop_levelset)[e->getID()] = levelsets;
+ }
+
+ auto mesh_prop_b = const_cast<MeshLib::Mesh&>(mesh)
+ .getProperties()
+ .template createNewPropertyVector<double>(
+ "aperture", MeshLib::MeshItemType::Cell);
+ mesh_prop_b->resize(mesh.getNumberOfElements());
+ auto mesh_prop_matid = mesh.getProperties().getPropertyVector<int>("MaterialIDs");
+ auto frac = _process_data._fracture_property.get();
+ for (MeshLib::Element const* e : _mesh.getElements())
+ {
+ if (e->getDimension() == DisplacementDim)
+ continue;
+ if ((*mesh_prop_matid)[e->getID()] != frac->mat_id)
+ continue;
+ ProcessLib::SpatialPosition x;
+ x.setElementID(e->getID());
+ (*mesh_prop_b)[e->getID()] = (*frac->aperture0)(0, x)[0];
+ }
+ _process_data._mesh_prop_b = mesh_prop_b;
+}
+
+
+template <int DisplacementDim>
+void SmallDeformationProcess<DisplacementDim>::postTimestepConcreteProcess(GlobalVector const& x)
+{
+ DBUG("PostTimestep SmallDeformationProcess.");
+
+ GlobalExecutor::executeMemberOnDereferenced(
+ &SmallDeformationLocalAssemblerInterface::postTimestep,
+ _local_assemblers, *_local_to_global_index_map, x);
+
+}
+
+
+// ------------------------------------------------------------------------------------
+// template instantiation
+// ------------------------------------------------------------------------------------
+template class SmallDeformationProcess<2>;
+
+} // namespace SmallDeformationWithLIE
+} // namespace ProcessLib
diff --git a/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.h b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.h
new file mode 100644
index 0000000000000000000000000000000000000000..33c408c7bc917b90bde06d3e14cd9d29d9724d31
--- /dev/null
+++ b/ProcessLib/SmallDeformationWithLIE/SmallDeformationProcess.h
@@ -0,0 +1,125 @@
+/**
+ * \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_SMALLDEFORMATION_WITH_LIE_PROCESS_H_
+#define PROCESSLIB_SMALLDEFORMATION_WITH_LIE_PROCESS_H_
+
+#include <cassert>
+
+#include "ProcessLib/Process.h"
+
+#include "ProcessLib/SmallDeformationWithLIE/LocalAssembler/CreateLocalAssemblers.h"
+#include "ProcessLib/SmallDeformationWithLIE/LocalAssembler/SmallDeformationLocalAssemblerInterface.h"
+
+#include "SmallDeformationProcessData.h"
+
+namespace ProcessLib
+{
+namespace SmallDeformationWithLIE
+{
+template <int DisplacementDim>
+class SmallDeformationProcess final : public Process
+{
+ using Base = Process;
+
+ static_assert(DisplacementDim==2,
+ "Currently SmallDeformationWithLIE::SmallDeformationProcess "
+ "supports only 2D.");
+
+public:
+ SmallDeformationProcess(
+ MeshLib::Mesh& mesh,
+ std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&&
+ jacobian_assembler,
+ std::vector<std::unique_ptr<ParameterBase>> const& parameters,
+ unsigned const integration_order,
+ std::vector<std::reference_wrapper<ProcessVariable>>&&
+ process_variables,
+ SmallDeformationProcessData<DisplacementDim>&& process_data,
+ SecondaryVariableCollection&& secondary_variables,
+ NumLib::NamedFunctionCaller&& named_function_caller);
+
+ //! \name ODESystem interface
+ //! @{
+
+ bool isLinear() const override { return false; }
+ //! @}
+
+private:
+ using LocalAssemblerInterface = SmallDeformationLocalAssemblerInterface;
+
+ void constructDofTable() override;
+
+ void initializeConcreteProcess(
+ NumLib::LocalToGlobalIndexMap const& dof_table,
+ MeshLib::Mesh const& mesh,
+ unsigned const integration_order) override;
+
+ void assembleConcreteProcess(const double t, GlobalVector const& x,
+ GlobalMatrix& M, GlobalMatrix& K,
+ GlobalVector& b) override
+ {
+ DBUG("Assemble SmallDeformationProcess.");
+
+ // Call global assembler for each local assembly item.
+ GlobalExecutor::executeMemberDereferenced(
+ _global_assembler, &VectorMatrixAssembler::assemble,
+ _local_assemblers, *_local_to_global_index_map, t, x, M, K, b);
+ }
+
+ void assembleWithJacobianConcreteProcess(
+ const double t, GlobalVector const& x, GlobalVector const& xdot,
+ const double dxdot_dx, const double dx_dx, GlobalMatrix& M,
+ GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) override
+ {
+ DBUG("AssembleWithJacobian SmallDeformationProcess.");
+
+ // Call global assembler for each local assembly item.
+ GlobalExecutor::executeMemberDereferenced(
+ _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
+ _local_assemblers, *_local_to_global_index_map, t, x, xdot,
+ dxdot_dx, dx_dx, M, K, b, Jac);
+ }
+
+ void preTimestepConcreteProcess(GlobalVector const& x, double const t,
+ double const dt) override
+ {
+ DBUG("PreTimestep SmallDeformationProcess.");
+
+ _process_data.dt = dt;
+ _process_data.t = t;
+
+ GlobalExecutor::executeMemberOnDereferenced(
+ &SmallDeformationLocalAssemblerInterface::preTimestep,
+ _local_assemblers, *_local_to_global_index_map, x, t, dt);
+ }
+
+ void postTimestepConcreteProcess(GlobalVector const& x) override;
+
+private:
+ SmallDeformationProcessData<DisplacementDim> _process_data;
+
+ std::vector<std::unique_ptr<LocalAssemblerInterface>> _local_assemblers;
+
+ std::unique_ptr<NumLib::LocalToGlobalIndexMap>
+ _local_to_global_index_map_single_component;
+
+ std::vector<MeshLib::Element*> _vec_matrix_elements;
+ std::vector<MeshLib::Element*> _vec_fracture_elements;
+ std::vector<MeshLib::Element*> _vec_fracture_matrix_elements;
+ std::vector<MeshLib::Node*> _vec_fracture_nodes;
+
+ std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_fracture_nodes;
+ std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_matrix_nodes;
+};
+
+} // namespace SmallDeformationWithLIE
+} // namespace ProcessLib
+
+#endif // PROCESSLIB_SMALLDEFORMATION_WITH_LIE_PROCESS_H_