From 9ac07ffc7454f0befd235d5d1c8e5a8badfe9f3a Mon Sep 17 00:00:00 2001
From: podlesny <podlesny@zedat.fu-berlin.de>
Date: Fri, 19 Feb 2021 16:30:21 +0100
Subject: [PATCH] save weights in programstate member
---
dune/tectonic/assemblers.cc | 3 ++-
dune/tectonic/assemblers.hh | 1 +
dune/tectonic/data-structures/program_state.hh | 7 ++++++-
3 files changed, 9 insertions(+), 2 deletions(-)
diff --git a/dune/tectonic/assemblers.cc b/dune/tectonic/assemblers.cc
index 37c8cb39..0f5ac990 100644
--- a/dune/tectonic/assemblers.cc
+++ b/dune/tectonic/assemblers.cc
@@ -144,6 +144,7 @@ template <class GridView, int dimension>
void MyAssembler<GridView, dimension>::assembleWeightedNormalStress(
BoundaryPatch<GridView> const &frictionalBoundary,
ScalarVector &weightedNormalStress,
+ ScalarVector &weights,
double youngModulus,
double poissonRatio,
Vector const &displacement) const {
@@ -159,7 +160,6 @@ void MyAssembler<GridView, dimension>::assembleWeightedNormalStress(
auto const nodalTractionAverage =
interpolateP0ToP1(frictionalBoundary, traction);
- ScalarVector weights;
{
NeumannBoundaryAssembler<Grid, typename ScalarVector::block_type>
frictionalBoundaryAssembler(
@@ -168,6 +168,7 @@ void MyAssembler<GridView, dimension>::assembleWeightedNormalStress(
vertexAssembler.assembleBoundaryFunctional(frictionalBoundaryAssembler,
weights, frictionalBoundary);
}
+
auto const normals = frictionalBoundary.getNormals();
for (size_t i = 0; i < vertexBasis.size(); ++i)
weightedNormalStress[i] =
diff --git a/dune/tectonic/assemblers.hh b/dune/tectonic/assemblers.hh
index e5354141..d7ca0d88 100644
--- a/dune/tectonic/assemblers.hh
+++ b/dune/tectonic/assemblers.hh
@@ -89,6 +89,7 @@ class MyAssembler {
void assembleWeightedNormalStress(
BoundaryPatch<GridView> const &frictionalBoundary,
ScalarVector &weightedNormalStress,
+ ScalarVector &weights,
double youngModulus,
double poissonRatio,
Vector const &displacement) const;
diff --git a/dune/tectonic/data-structures/program_state.hh b/dune/tectonic/data-structures/program_state.hh
index 4bc169ec..d4573f9b 100644
--- a/dune/tectonic/data-structures/program_state.hh
+++ b/dune/tectonic/data-structures/program_state.hh
@@ -74,6 +74,7 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
v(bodyCount_),
a(bodyCount_),
alpha(bodyCount_),
+ weights(bodyCount_),
weightedNormalStress(bodyCount_) {
for (size_t i=0; i<bodyCount_; i++) {
size_t leafVertexCount = leafVertexCounts[i];
@@ -82,6 +83,7 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
v[i].resize(leafVertexCount),
a[i].resize(leafVertexCount),
alpha[i].resize(leafVertexCount),
+ weights[i].resize(leafVertexCount),
weightedNormalStress[i].resize(leafVertexCount),
bodyStates[i] = new BodyState(&u[i], &v[i], &a[i], &alpha[i], &weightedNormalStress[i]);
@@ -205,12 +207,14 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
const auto& body = contactNetwork.body(nonmortarIdx);
ScalarVector frictionBoundaryStress(weightedNormalStress[nonmortarIdx].size());
+ ScalarVector nmWeights(frictionBoundaryStress);
body->assembler()->assembleWeightedNormalStress(
- contactCoupling->nonmortarBoundary(), frictionBoundaryStress, body->data()->getYoungModulus(),
+ contactCoupling->nonmortarBoundary(), frictionBoundaryStress, nmWeights, body->data()->getYoungModulus(),
body->data()->getPoissonRatio(), u[nonmortarIdx]);
weightedNormalStress[nonmortarIdx] += frictionBoundaryStress;
+ weights[nonmortarIdx] += nmWeights;
}
std::cout << "solving linear problem for a..." << std::endl;
@@ -231,6 +235,7 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
std::vector<Vector> v;
std::vector<Vector> a;
std::vector<ScalarVector> alpha;
+ std::vector<ScalarVector> weights;
std::vector<ScalarVector> weightedNormalStress;
double relativeTime;
double relativeTau;
--
GitLab