From 0fcefc8100d4ffe5c79b54e871cb2d710f857730 Mon Sep 17 00:00:00 2001
From: podlesny <podlesny@mi.fu-berlin.de>
Date: Thu, 25 Jul 2019 14:32:43 +0200
Subject: [PATCH] .
---
src/data-structures/contactnetwork.cc | 20 +++++++++++-------
src/multi-body-problem-2D.cfg | 2 +-
src/multi-body-problem-data/bc.hh | 4 ++--
src/multi-body-problem.cfg | 2 +-
src/solverfactorytest.cc | 21 ++++++++++++-------
src/time-stepping/rate/newmark.cc | 6 +++---
.../state/ageinglawstateupdater.cc | 2 +-
7 files changed, 33 insertions(+), 24 deletions(-)
diff --git a/src/data-structures/contactnetwork.cc b/src/data-structures/contactnetwork.cc
index a25f7e49..84517f20 100644
--- a/src/data-structures/contactnetwork.cc
+++ b/src/data-structures/contactnetwork.cc
@@ -26,7 +26,7 @@ ContactNetwork<HostGridType, VectorType>::ContactNetwork(
couplings_(nCouplings),
frictionBoundaries_(nBodies),
stateEnergyNorms_(nBodies),
- nBodyAssembler_(nBodies, nCouplings) //, false, 0.0)
+ nBodyAssembler_(nBodies, nCouplings, false, 0.0)
{}
template <class HostGridType, class VectorType>
@@ -324,13 +324,17 @@ void ContactNetwork<HostGridType, VectorType>::totalNodes(
std::vector<std::shared_ptr<typename LeafBody::BoundaryCondition>> boundaryConditions;
body->boundaryConditions(tag, boundaryConditions);
- const int idxBackup = idx;
- for (size_t bc = 0; bc<boundaryConditions.size(); bc++) {
- const auto& nodes = boundaryConditions[bc]->boundaryNodes();
- for (size_t j=0; j<nodes->size(); j++, idx++)
- for (int k=0; k<dim; k++)
- totalNodes[idx][k] = (*nodes)[j][k];
- idx = (bc==boundaryConditions.size()-1 ? idx : idxBackup);
+ if (boundaryConditions.size()>0) {
+ const int idxBackup = idx;
+ for (size_t bc = 0; bc<boundaryConditions.size(); bc++) {
+ const auto& nodes = boundaryConditions[bc]->boundaryNodes();
+ for (size_t j=0; j<nodes->size(); j++, idx++)
+ for (int k=0; k<dim; k++)
+ totalNodes[idx][k] = (*nodes)[j][k];
+ idx = (bc==boundaryConditions.size()-1 ? idx : idxBackup);
+ }
+ } else {
+ idx += body->nVertices();
}
}
}
diff --git a/src/multi-body-problem-2D.cfg b/src/multi-body-problem-2D.cfg
index a4a77808..a308c80b 100644
--- a/src/multi-body-problem-2D.cfg
+++ b/src/multi-body-problem-2D.cfg
@@ -1,6 +1,6 @@
# -*- mode:conf -*-
[boundary.friction]
-smallestDiameter = 0.25 # 2e-3 [m]
+smallestDiameter = 0.01 # 2e-3 [m]
[timeSteps]
refinementTolerance = 1e-5 # 1e-5
diff --git a/src/multi-body-problem-data/bc.hh b/src/multi-body-problem-data/bc.hh
index c7c919ae..b9b7ecda 100644
--- a/src/multi-body-problem-data/bc.hh
+++ b/src/multi-body-problem-data/bc.hh
@@ -9,12 +9,12 @@ class VelocityDirichletCondition
// Assumed to vanish at time zero
double const finalVelocity = -5e-5;
- std::cout << "VelocityDirichletCondition::evaluate()" << std::endl;
+ /*std::cout << "VelocityDirichletCondition::evaluate()" << std::endl;
if (relativeTime <= 0.1)
std::cout << "- loading phase" << std::endl;
else
- std::cout << "- final velocity reached" << std::endl;
+ std::cout << "- final velocity reached" << std::endl;*/
y = (relativeTime <= 0.1)
? finalVelocity * (1.0 - std::cos(relativeTime * M_PI / 0.1)) / 2.0
diff --git a/src/multi-body-problem.cfg b/src/multi-body-problem.cfg
index 7188eed1..e8a3e97b 100644
--- a/src/multi-body-problem.cfg
+++ b/src/multi-body-problem.cfg
@@ -11,7 +11,7 @@ vtk.write = true
[problem]
finalTime = 100 # [s] #1000
-bodyCount = 2
+bodyCount = 4
[body]
bulkModulus = 0.5e5 # [Pa]
diff --git a/src/solverfactorytest.cc b/src/solverfactorytest.cc
index 19b68965..1d0f983b 100644
--- a/src/solverfactorytest.cc
+++ b/src/solverfactorytest.cc
@@ -57,7 +57,7 @@
const int dim = MY_DIM;
-const int timeSteps = 1;
+const int timeSteps = 5;
Dune::ParameterTree parset;
size_t bodyCount;
@@ -86,6 +86,11 @@ void solveProblem(const ContactNetwork& contactNetwork,
//using LocalSolver = LocalBisectionSolver;
//using Linearization = Linearization<Functional, BitVector>;
+ /*print(ignore, "ignore:");
+ for (size_t i=0; i<x.size(); i++) {
+ std::cout << x[i] << std::endl;
+ }*/
+
// set up reference solver
Vector refX = x;
using ContactFunctional = Dune::TNNMG::BoxConstrainedQuadraticFunctional<Matrix&, Vector&, Vector&, Vector&, double>;
@@ -438,8 +443,8 @@ void run(Updaters& updaters, ContactNetwork& contactNetwork,
nBodyAssembler.concatenateVectors(u0, totalu0);
nBodyAssembler.concatenateVectors(v0, totalv0);
- updaters.rate_->velocityObstacles(totalu0, lower, totalv0, vLower);
- updaters.rate_->velocityObstacles(totalu0, upper, totalv0, vUpper);
+ //updaters.rate_->velocityObstacles(totalu0, lower, totalv0, vLower);
+ //updaters.rate_->velocityObstacles(totalu0, upper, totalv0, vUpper);
const auto& errorNorms = contactNetwork.stateEnergyNorms();
@@ -465,7 +470,7 @@ void run(Updaters& updaters, ContactNetwork& contactNetwork,
Vector totalVelocityIterate;
nBodyAssembler.nodalToTransformed(velocityIterates, totalVelocityIterate);
- solveProblem(contactNetwork, bilinearForm, totalRhs, totalVelocityIterate, totalDirichletNodes, vLower, vUpper);
+ solveProblem(contactNetwork, bilinearForm, totalRhs, totalVelocityIterate, totalDirichletNodes, lower, upper);
nBodyAssembler.postprocess(totalVelocityIterate, velocityIterates);
@@ -496,8 +501,8 @@ void run(Updaters& updaters, ContactNetwork& contactNetwork,
if (alpha[i].size()==0 || newAlpha[i].size()==0)
continue;
- print(alpha[i], "alpha i:");
- print(newAlpha[i], "new alpha i:");
+ //print(alpha[i], "alpha i:");
+ //print(newAlpha[i], "new alpha i:");
if (errorNorms[i]->diff(alpha[i], newAlpha[i]) >= fixedPointTolerance) {
breakCriterion = false;
std::cout << "fixedPoint error: " << errorNorms[i]->diff(alpha[i], newAlpha[i]) << std::endl;
@@ -674,10 +679,10 @@ int main(int argc, char *argv[]) { try {
// Set up TNNMG solver
// -------------------
- BitVector totalDirichletNodes;
+ /*BitVector totalDirichletNodes;
contactNetwork.totalNodes("dirichlet", totalDirichletNodes);
- //print(totalDirichletNodes, "totalDirichletNodes:");
+ print(totalDirichletNodes, "totalDirichletNodes:");*/
//using Functional = Functional<Matrix&, Vector&, ZeroNonlinearity&, Vector&, Vector&, field_type>;
//using Functional = Functional<Matrix&, Vector&, GlobalFriction<Matrix, Vector>&, Vector&, Vector&, field_type>;
diff --git a/src/time-stepping/rate/newmark.cc b/src/time-stepping/rate/newmark.cc
index a3b36385..0ec59012 100644
--- a/src/time-stepping/rate/newmark.cc
+++ b/src/time-stepping/rate/newmark.cc
@@ -95,20 +95,20 @@ void Newmark<Vector, Matrix, BoundaryFunctions, BoundaryNodes>::setup(
double dirichletValue;
bodyDirichletFunctions[bc]->evaluate(relativeTime, dirichletValue);
- std::cout << "dirichletValue: " << dirichletValue << std::endl;
+ //std::cout << "dirichletValue: " << dirichletValue << std::endl;
for (size_t k=0; k<bcDirichletNodes.size() ; ++k) {
for (size_t j=0; j<dim; ++j) {
if (bcDirichletNodes[k][j]) {
- std::cout << k << " " << j << std::endl;
+ //std::cout << k << " " << j << std::endl;
bodyIterate[k][j] = dirichletValue;
}
}
}
}
}
- print(iterate, "iterate:");
+ //sprint(iterate, "iterate:");
}
template <class Vector, class Matrix, class BoundaryFunctions, class BoundaryNodes>
diff --git a/src/time-stepping/state/ageinglawstateupdater.cc b/src/time-stepping/state/ageinglawstateupdater.cc
index 796f8878..1a7dcea3 100644
--- a/src/time-stepping/state/ageinglawstateupdater.cc
+++ b/src/time-stepping/state/ageinglawstateupdater.cc
@@ -66,7 +66,7 @@ template <class ScalarVector, class Vector>
void AgeingLawStateUpdater<ScalarVector, Vector>::extractAlpha(
ScalarVector& alpha) {
- std::cout << "alpha size: " << alpha.size() << " nodes_.size() " << nodes_.size() << std::endl;
+ //std::cout << "alpha size: " << alpha.size() << " nodes_.size() " << nodes_.size() << std::endl;
if (alpha.size() != nodes_.size()) {
alpha.resize(nodes_.size());
}
--
GitLab