From 50a9b11676f5dd23822227cd92b4476dd2235f19 Mon Sep 17 00:00:00 2001
From: "Praetorius, Simon" <simon.praetorius@tu-dresden.de>
Date: Wed, 28 Aug 2019 11:10:56 +0200
Subject: [PATCH] added stress tensor part for non-constant viscosity
---
src/amdis/localoperators/StokesOperator.hpp | 22 ++++++++++++++++++++-
1 file changed, 21 insertions(+), 1 deletion(-)
diff --git a/src/amdis/localoperators/StokesOperator.hpp b/src/amdis/localoperators/StokesOperator.hpp
index a8f03f03..1591f37c 100644
--- a/src/amdis/localoperators/StokesOperator.hpp
+++ b/src/amdis/localoperators/StokesOperator.hpp
@@ -31,8 +31,9 @@ namespace AMDiS
static_assert( Size_v<typename ViscosityExpr::Range> == 1, "Viscosity must be of scalar type." );
public:
- GridFunctionOperator(tag::stokes, ViscosityExpr const& expr)
+ GridFunctionOperator(tag::stokes, ViscosityExpr const& expr, bool constViscosity = false)
: Super(expr, 1)
+ , constViscosity_(constViscosity)
{}
template <class CG, class Node, class Mat>
@@ -100,6 +101,22 @@ namespace AMDiS
}
}
+ if (!constViscosity_) {
+ // <viscosity * d_i u_j, d_j v_i>
+ for (std::size_t i = 0; i < numVelocityLocalFE; ++i) {
+ for (std::size_t kj = 0; kj < CG::dow; ++kj) {
+ const auto value_i_kj = vel_factor * gradients[i][kj];
+ for (std::size_t j = 0; j < numVelocityLocalFE; ++j) {
+ for (std::size_t ki = 0; ki < CG::dow; ++ki) {
+ const auto local_ki = tree.child(_0,ki).localIndex(i);
+ const auto local_kj = tree.child(_0,kj).localIndex(j);
+ elementMatrix[local_ki][local_kj] += value_i_kj * gradients[j][ki];
+ }
+ }
+ }
+ }
+ }
+
// <p, div(v)> + <div(u), q>
for (std::size_t i = 0; i < numVelocityLocalFE; ++i) {
for (std::size_t j = 0; j < numPressureLocalFE; ++j) {
@@ -116,6 +133,9 @@ namespace AMDiS
}
} // end calculateElementMatrix
+
+ private:
+ bool constViscosity_;
};
/** @} **/
--
GitLab