diff --git a/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_solid_density.md b/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_solid_density.md
new file mode 100644
index 0000000000000000000000000000000000000000..4cc76c579164f400870a8eb0ac4d6c8d58e2c13e
--- /dev/null
+++ b/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_solid_density.md
@@ -0,0 +1 @@
+\copydoc ProcessLib::SmallDeformation::SmallDeformationProcessData::solid_density
diff --git a/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_specific_body_force.md b/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_specific_body_force.md
new file mode 100644
index 0000000000000000000000000000000000000000..a481e852d81314336ca7ec474924e2763c2209d3
--- /dev/null
+++ b/Documentation/ProjectFile/prj/processes/process/SMALL_DEFORMATION/t_specific_body_force.md
@@ -0,0 +1 @@
+\copydoc ProcessLib::SmallDeformation::SmallDeformationProcessData::specific_body_force
diff --git a/ProcessLib/SmallDeformation/CreateSmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/CreateSmallDeformationProcess.cpp
index ceb84f4daed87b3042c703bbadf2a4b382acdd17..8db7dc6dfc234d1d6b49bb569729339b8da57f20 100644
--- a/ProcessLib/SmallDeformation/CreateSmallDeformationProcess.cpp
+++ b/ProcessLib/SmallDeformation/CreateSmallDeformationProcess.cpp
@@ -97,8 +97,32 @@ createSmallDeformationProcess(
type.c_str());
}
+ // Solid density
+ auto& solid_density = findParameter<double>(
+ config,
+ //! \ogs_file_param_special{prj__processes__process__SMALL_DEFORMATION__solid_density}
+ "solid_density", parameters, 1);
+ DBUG("Use \'%s\' as solid density parameter.", solid_density.name.c_str());
+
+ // Specific body force
+ Eigen::Matrix<double, DisplacementDim, 1> specific_body_force;
+ {
+ std::vector<double> const b =
+ //! \ogs_file_param{prj__processes__process__SMALL_DEFORMATION__specific_body_force}
+ config.getConfigParameter<std::vector<double>>(
+ "specific_body_force");
+ if (b.size() != DisplacementDim)
+ OGS_FATAL(
+ "The size of the specific body force vector does not match the "
+ "displacement dimension. Vector size is %d, displacement "
+ "dimension is %d",
+ b.size(), DisplacementDim);
+
+ std::copy_n(b.data(), b.size(), specific_body_force.data());
+ }
+
SmallDeformationProcessData<DisplacementDim> process_data{
- std::move(material)};
+ std::move(material), solid_density, specific_body_force};
SecondaryVariableCollection secondary_variables;
diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index c2865ed1ebd67716a4302f41ae291cf5d41cf47c..7d20d5ee34a2cb13ab8dcf17b929cfbcd451f7e4 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
@@ -185,6 +185,18 @@ public:
auto const& N = _ip_data[ip].N;
auto const& dNdx = _ip_data[ip].dNdx;
+ typename ShapeMatricesType::template MatrixType<DisplacementDim,
+ displacement_size>
+ N_u_op = ShapeMatricesType::template MatrixType<
+ DisplacementDim,
+ displacement_size>::Zero(DisplacementDim,
+ displacement_size);
+ for (int i = 0; i < DisplacementDim; ++i)
+ N_u_op
+ .template block<1, displacement_size / DisplacementDim>(
+ i, i * displacement_size / DisplacementDim)
+ .noalias() = N;
+
auto const x_coord =
interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
_element, N);
@@ -215,7 +227,10 @@ public:
KelvinMatrixType<DisplacementDim> C;
std::tie(sigma, state, C) = std::move(*solution);
- local_b.noalias() -= B.transpose() * sigma * w;
+ auto const rho = _process_data.solid_density(t, x_position)[0];
+ auto const& b = _process_data.specific_body_force;
+ local_b.noalias() -=
+ (B.transpose() * sigma - N_u_op.transpose() * rho * b) * w;
local_Jac.noalias() += B.transpose() * C * B * w;
}
}
@@ -493,6 +508,9 @@ private:
MeshLib::Element const& _element;
SecondaryData<typename ShapeMatrices::ShapeType> _secondary_data;
bool const _is_axially_symmetric;
+
+ static const int displacement_size =
+ ShapeFunction::NPOINTS * DisplacementDim;
};
} // namespace SmallDeformation
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcessData.h b/ProcessLib/SmallDeformation/SmallDeformationProcessData.h
index 3aca5d5630dc0f8263ec91a3bdf41523a00dd6c2..115b0f286cc1b154996f053032b1d7e3fc435d51 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcessData.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcessData.h
@@ -26,13 +26,22 @@ struct SmallDeformationProcessData
{
SmallDeformationProcessData(
std::unique_ptr<MaterialLib::Solids::MechanicsBase<DisplacementDim>>&&
- material)
- : material{std::move(material)}
+ material,
+ Parameter<double> const& solid_density_,
+ Eigen::Matrix<double, DisplacementDim, 1>
+ specific_body_force_)
+ : material{std::move(material)},
+ solid_density(solid_density_),
+ specific_body_force(std::move(specific_body_force_))
{
}
SmallDeformationProcessData(SmallDeformationProcessData&& other)
- : material{std::move(other.material)}, dt{other.dt}, t{other.t}
+ : material{std::move(other.material)},
+ solid_density(other.solid_density),
+ specific_body_force(other.specific_body_force),
+ dt{other.dt},
+ t{other.t}
{
}
@@ -47,6 +56,12 @@ struct SmallDeformationProcessData
std::unique_ptr<MaterialLib::Solids::MechanicsBase<DisplacementDim>>
material;
+ /// Solid's density. A scalar quantity, Parameter<double>.
+ Parameter<double> const& solid_density;
+ /// Specific body forces applied to the solid.
+ /// It is usually used to apply gravitational forces.
+ /// A vector of displacement dimension's length.
+ Eigen::Matrix<double, DisplacementDim, 1> const specific_body_force;
double dt = 0;
double t = 0;
};
diff --git a/ProcessLib/SmallDeformation/Tests.cmake b/ProcessLib/SmallDeformation/Tests.cmake
index f7c32f4133dd373beb1bd8d989a2cca45f0065d6..8d639fc02fb042af27fc0fabe9e5136c5147529d 100644
--- a/ProcessLib/SmallDeformation/Tests.cmake
+++ b/ProcessLib/SmallDeformation/Tests.cmake
@@ -327,3 +327,15 @@ AddTest(
ref_axisymmetric_sphere_pcs_0_ts_100_t_1.000000.vtu axisymmetric_sphere_pcs_0_ts_100_t_1.000000.vtu sigma_yy sigma_yy
ref_axisymmetric_sphere_pcs_0_ts_100_t_1.000000.vtu axisymmetric_sphere_pcs_0_ts_100_t_1.000000.vtu sigma_zz sigma_zz
)
+
+# Two materials in gravity field
+AddTest(
+ NAME Parallel_Mechanics_SDL_two_material_gravity
+ PATH Mechanics/Linear
+ EXECUTABLE_ARGS two_material_gravity.prj
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ ABSTOL 1e-16 RELTOL 1e-16
+ DIFF_DATA
+ expected_two_material_gravity.vtu two_material_gravity_pcs_0_ts_1_t_1.000000.vtu displacement displacement
+)
diff --git a/Tests/Data b/Tests/Data
index c8b8a5628f9e0de155ff681780e485d8ba0dd18b..e31b2e066e2705a3857c1ffbc1cd3fffbdd1d379 160000
--- a/Tests/Data
+++ b/Tests/Data
@@ -1 +1 @@
-Subproject commit c8b8a5628f9e0de155ff681780e485d8ba0dd18b
+Subproject commit e31b2e066e2705a3857c1ffbc1cd3fffbdd1d379