From c38d1abd8e671cb71f54dd7f4fcd2633fea128b2 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <dmitri.naumov@ufz.de>
Date: Fri, 8 Sep 2017 18:15:38 +0200
Subject: [PATCH] [PL] SD: Add body forces to the FEM.
---
.../SmallDeformation/SmallDeformationFEM.h | 20 ++++++++++++++++++-
1 file changed, 19 insertions(+), 1 deletion(-)
diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index c2865ed1ebd..7d20d5ee34a 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
--
GitLab