diff --git a/dune/tectonic/globalratestatefriction.hh b/dune/tectonic/globalratestatefriction.hh
index edfb05e056b9a88ebf265bb4f123dfab04c85b1f..078bd85dc9cf1ab19e5a7427ecd160a91dc3f9aa 100644
--- a/dune/tectonic/globalratestatefriction.hh
+++ b/dune/tectonic/globalratestatefriction.hh
@@ -88,7 +88,7 @@ class GlobalRateStateFriction : public GlobalFriction<Matrix, Vector> {
localToGlobal_.emplace_back(i);
restrictions_.emplace_back(weights[bodyIdx][i], weightedNormalStress[bodyIdx][i],
- couplings[i]->frictionData()(geoToPoint(it->geometry())));
+ couplings[j]->frictionData()(geoToPoint(it->geometry())));
break;
}
diff --git a/src/data-structures/levelcontactnetwork.cc b/src/data-structures/levelcontactnetwork.cc
index 9709a867724b7b30c9bcfd7de7815511aafae33b..aa3365139c3820c8f37b07887c070a904fc25b52 100644
--- a/src/data-structures/levelcontactnetwork.cc
+++ b/src/data-structures/levelcontactnetwork.cc
@@ -27,6 +27,8 @@ LevelContactNetwork<GridType, dimension>::LevelContactNetwork(
leafViews_(nBodies),
levelViews_(nBodies),
leafVertexCounts_(nBodies),
+ frictionBoundaries_(nCouplings),
+ stateEnergyNorms_(nBodies),
nBodyAssembler_(nBodies, nCouplings)
{}
@@ -34,6 +36,7 @@ template <class GridType, int dimension>
LevelContactNetwork<GridType, dimension>::~LevelContactNetwork() {
for (size_t i=0; i<nBodies(); i++) {
delete frictionBoundaries_[i];
+ delete frictionBoundaryMass_[i];
delete stateEnergyNorms_[i];
}
}
@@ -46,33 +49,37 @@ void LevelContactNetwork<GridType, dimension>::assemble() {
bodies_[i]->assemble();
//printBasisDofLocation(bodies_[i]->assembler()->vertexBasis); //TODO remove after debugging
- frictionBoundaries_[i] = new BoundaryPatch(bodies_[i]->leafView(),false);
+ frictionBoundaries_[i] = new BoundaryPatch(bodies_[i]->leafView(), false);
}
// set up dune-contact nBodyAssembler
nBodyAssembler_.setGrids(deformedGrids_);
for (size_t i=0; i<nCouplings(); i++) {
- auto& coupling = couplings_[i];
- nBodyAssembler_.setCoupling(*coupling, i);
+ nBodyAssembler_.setCoupling(*couplings_[i], i);
+ }
+
+ nBodyAssembler_.assembleTransferOperator();
+ nBodyAssembler_.assembleObstacle();
+ for (size_t i=0; i<nCouplings(); i++) {
+ auto& coupling = couplings_[i];
const auto& contactCoupling = nBodyAssembler_.getContactCouplings()[i]; // dual mortar object holding boundary patches
const auto nonmortarIdx = coupling->gridIdx_[0];
- const auto mortarIdx = coupling->gridIdx_[1];
+ //const auto mortarIdx = coupling->gridIdx_[1];
frictionBoundaries_[nonmortarIdx]->addPatch(contactCoupling->nonmortarBoundary());
- frictionBoundaries_[mortarIdx]->addPatch(contactCoupling->mortarBoundary());
+ //frictionBoundaries_[mortarIdx]->addPatch(contactCoupling->mortarBoundary());
}
- nBodyAssembler_.assembleTransferOperator();
- nBodyAssembler_.assembleObstacle();
-
// assemble state energy norm
+ frictionBoundaryMass_.resize(nBodies());
for (size_t i=0; i<nBodies(); i++) {
- ScalarMatrix relativeFrictionBoundaryMass;
- bodies_[i]->assembler()->assembleFrictionalBoundaryMass(*frictionBoundaries_[i], relativeFrictionBoundaryMass);
- relativeFrictionBoundaryMass /= frictionBoundaries_[i]->area(); // TODO: weight by individual friction patches?
- stateEnergyNorms_[i] = new EnergyNorm<ScalarMatrix, ScalarVector>(relativeFrictionBoundaryMass);
+ frictionBoundaryMass_[i] = new ScalarMatrix();
+ bodies_[i]->assembler()->assembleFrictionalBoundaryMass(*frictionBoundaries_[i], *frictionBoundaryMass_[i]);
+ *frictionBoundaryMass_[i] /= frictionBoundaries_[i]->area(); // TODO: weight by individual friction patches?
+
+ stateEnergyNorms_[i] = new EnergyNorm<ScalarMatrix, ScalarVector>(*frictionBoundaryMass_[i]);
}
}
@@ -103,12 +110,15 @@ void LevelContactNetwork<GridType, dimension>::assembleFriction(
globalFriction_ = std::make_shared<GlobalRateStateFriction<
Matrix, Vector, TruncatedRateState, DeformedGridType>>(
nBodyAssembler_.getContactCouplings(), couplings_, weights, weightedNormalStress);
+ break;
case Config::Regularised:
globalFriction_ = std::make_shared<GlobalRateStateFriction<
Matrix, Vector, RegularisedRateState, DeformedGridType>>(
nBodyAssembler_.getContactCouplings(), couplings_, weights, weightedNormalStress);
+ break;
default:
assert(false);
+ break;
}
}
@@ -266,6 +276,7 @@ void LevelContactNetwork<GridType, dimension>::frictionNodes(std::vector<const D
template <class GridType, int dimension>
auto LevelContactNetwork<GridType, dimension>::stateEnergyNorms() const
-> const StateEnergyNorms& {
+
return stateEnergyNorms_;
}
diff --git a/src/data-structures/levelcontactnetwork.hh b/src/data-structures/levelcontactnetwork.hh
index 4721208b1e6b744f59cccf96e8c05f145ed6e4c3..97792f2ee6dcd3ff2483e81e1eba0d053a1b3322 100644
--- a/src/data-structures/levelcontactnetwork.hh
+++ b/src/data-structures/levelcontactnetwork.hh
@@ -143,6 +143,7 @@ template <class GridType, int dimension> class LevelContactNetwork {
std::vector<size_t> leafVertexCounts_;
std::vector<DeformedLeafBoundaryPatch*> frictionBoundaries_; // union of all boundary patches per body
+ std::vector<ScalarMatrix*> frictionBoundaryMass_;
StateEnergyNorms stateEnergyNorms_;
NBodyAssembler nBodyAssembler_;
diff --git a/src/data-structures/program_state.hh b/src/data-structures/program_state.hh
index ec0311fef46e1f2fdf9cf03505c69e79c5b326b1..298f40f7e14fbc10b17217fd8fd59da4be765bb1 100644
--- a/src/data-structures/program_state.hh
+++ b/src/data-structures/program_state.hh
@@ -143,19 +143,19 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
}
// print problem
- print(bilinearForm, "bilinearForm");
+ /* print(bilinearForm, "bilinearForm");
print(totalRhs, "totalRhs");
print(_dirichletNodes, "ignore");
print(totalObstacles, "totalObstacles");
print(lower, "lower");
- print(upper, "upper");
+ print(upper, "upper");*/
// set up functional and TNMMG solver
using Functional = Functional<Matrix&, Vector&, ZeroNonlinearity&, Vector&, Vector&, typename Matrix::field_type>;
Functional J(bilinearForm, totalRhs, ZeroNonlinearity(), lower, upper);
using Factory = SolverFactory<Functional, Dune::BitSetVector<dims>>;
- Factory factory(parset.sub("solver.tnnmg"), J);
+ Factory factory(parset.sub("solver.tnnmg"), J, _dirichletNodes);
/* std::vector<Dune::BitSetVector<dims>> bodyDirichletNodes;
nBodyAssembler.postprocess(_dirichletNodes, bodyDirichletNodes);
@@ -168,7 +168,6 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
print(totalRhs, "totalRhs: ");*/
auto tnnmgStep = factory.step();
- tnnmgStep->setIgnore(_dirichletNodes);
factory.setProblem(totalX);
const EnergyNorm<Matrix, Vector> norm(bilinearForm);
@@ -209,7 +208,7 @@ template <class VectorTEMPLATE, class ScalarVectorTEMPLATE> class ProgramState {
bilinearForm.mmv(oldStep->getSol(), oldRes);
std::cout << "Old Res - energy norm: " << norm.operator ()(oldRes) << std::endl;*/
- print(tnnmgStep->getSol(), "TNNMG Solution: ");
+ // print(tnnmgStep->getSol(), "TNNMG Solution: ");
/* print(oldStep->getSol(), "Old Solution: ");
auto diff = tnnmgStep->getSol();
diff -= oldStep->getSol();
diff --git a/src/multi-body-problem.cc b/src/multi-body-problem.cc
index 709adff5196f71f5c17d02480a8d12f19f60122f..fda8fc8f1308835c69386b6650a0dedc7ac13dea 100644
--- a/src/multi-body-problem.cc
+++ b/src/multi-body-problem.cc
@@ -215,8 +215,8 @@ int main(int argc, char *argv[]) {
vertexCoords[vertexMapper.index(v)] = geoToPoint(v.geometry());
}
- typename LevelContactNetwork::BoundaryPatches frictionBoundaries;
- levelContactNetwork.boundaryPatches("friction", frictionBoundaries);
+ //typename LevelContactNetwork::BoundaryPatches frictionBoundaries;
+ //levelContactNetwork.boundaryPatches("friction", frictionBoundaries);
/*
auto dataWriter =
@@ -273,7 +273,7 @@ int main(int argc, char *argv[]) {
// -------------------
auto& nBodyAssembler = levelContactNetwork.nBodyAssembler();
- Dune::BitSetVector<dims> totalDirichletNodes;
+ BitVector totalDirichletNodes;
levelContactNetwork.totalNodes("dirichlet", totalDirichletNodes);
using Functional = Functional<Matrix&, Vector&, GlobalFriction<Matrix, Vector>&, Vector&, Vector&, field_type>;
@@ -338,7 +338,7 @@ int main(int argc, char *argv[]) {
levelContactNetwork.externalForces(externalForces);
AdaptiveTimeStepper<NonlinearFactory, std::decay_t<decltype(nBodyAssembler)>, Updaters, EnergyNorm<ScalarMatrix, ScalarVector>>
- adaptiveTimeStepper(parset, nBodyAssembler, globalFriction, frictionNodes, current,
+ adaptiveTimeStepper(parset, nBodyAssembler, totalDirichletNodes, globalFriction, frictionNodes, current,
programState.relativeTime, programState.relativeTau,
externalForces, stateEnergyNorms, mustRefine);
diff --git a/src/spatial-solving/fixedpointiterator.cc b/src/spatial-solving/fixedpointiterator.cc
index 5ca8492670561d87d71ec486ea64122642cbcb30..ccf52350a9303cf8792b011160ffe59bba249cd0 100644
--- a/src/spatial-solving/fixedpointiterator.cc
+++ b/src/spatial-solving/fixedpointiterator.cc
@@ -27,8 +27,10 @@
#include "fixedpointiterator.hh"
+#include "../utils/tobool.hh"
#include "../utils/debugutils.hh"
+
void FixedPointIterationCounter::operator+=(
FixedPointIterationCounter const &other) {
iterations += other.iterations;
@@ -39,11 +41,13 @@ template <class Factory, class NBodyAssembler, class Updaters, class ErrorNorm>
FixedPointIterator<Factory, NBodyAssembler, Updaters, ErrorNorm>::FixedPointIterator(
Dune::ParameterTree const &parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
const std::vector<const ErrorNorm* >& errorNorms)
: parset_(parset),
nBodyAssembler_(nBodyAssembler),
+ ignoreNodes_(ignoreNodes),
globalFriction_(globalFriction),
bodywiseNonmortarBoundaries_(bodywiseNonmortarBoundaries),
fixedPointMaxIterations_(parset.get<size_t>("v.fpi.maximumIterations")),
@@ -94,7 +98,7 @@ FixedPointIterator<Factory, NBodyAssembler, Updaters, ErrorNorm>::run(
// set up functional and TNMMG solver
Functional J(bilinearForm, totalRhs, globalFriction_, lower, upper);
- Factory solverFactory(parset_.sub("solver.tnnmg"), J);
+ Factory solverFactory(parset_.sub("solver.tnnmg"), J, ignoreNodes_);
auto step = solverFactory.step();
EnergyNorm<Matrix, Vector> energyNorm(bilinearForm);
@@ -183,10 +187,24 @@ void FixedPointIterator<Factory, NBodyAssembler, Updaters, ErrorNorm>::relativeV
const size_t nBodies = nBodyAssembler_.nGrids();
// const auto& contactCouplings = nBodyAssembler_.getContactCouplings();
-
+ size_t globalIdx = 0;
v_rel.resize(nBodies);
+ for (size_t bodyIdx=0; bodyIdx<nBodies; bodyIdx++) {
+ const auto& nonmortarBoundary = *bodywiseNonmortarBoundaries_[bodyIdx];
+ auto& v_rel_ = v_rel[bodyIdx];
+
+ v_rel_.resize(nonmortarBoundary.size());
+
+ for (size_t i=0; i<v_rel_.size(); i++) {
+ if (toBool(nonmortarBoundary[i])) {
+ v_rel_[i] = v[globalIdx];
+ }
+ globalIdx++;
+ }
+ }
+
- /* for (size_t bodyIdx=0; bodyIdx<nBodies; bodyIdx++) {
+ /*
boundaryNodes =
const auto gridView = contactCouplings[couplingIndices[0]]->nonmortarBoundary().gridView();
diff --git a/src/spatial-solving/fixedpointiterator.hh b/src/spatial-solving/fixedpointiterator.hh
index 078b1d2051c539538b4cc3b0eebfb0f42851f6fe..fa93c0870555ca106090401e3c6f1ec582d999bc 100644
--- a/src/spatial-solving/fixedpointiterator.hh
+++ b/src/spatial-solving/fixedpointiterator.hh
@@ -33,6 +33,7 @@ class FixedPointIterator {
const static int dims = Vector::block_type::dimension;
+ using IgnoreVector = typename Factory::BitVector;
// using Nonlinearity = typename ConvexProblem::NonlinearityType;
// using DeformedGrid = typename Factory::DeformedGrid;
@@ -46,6 +47,7 @@ class FixedPointIterator {
public:
FixedPointIterator(const Dune::ParameterTree& parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
const std::vector<const ErrorNorm* >& errorNorms);
@@ -58,6 +60,8 @@ class FixedPointIterator {
private:
const Dune::ParameterTree& parset_;
NBodyAssembler& nBodyAssembler_;
+ const IgnoreVector& ignoreNodes_;
+
GlobalFriction& globalFriction_;
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries_;
diff --git a/src/spatial-solving/solverfactory.cc b/src/spatial-solving/solverfactory.cc
index 25bc75755807abb46b4201ed8f25aa319732cc21..4804b78837521099d09ee4cab006c045e77b9266 100644
--- a/src/spatial-solving/solverfactory.cc
+++ b/src/spatial-solving/solverfactory.cc
@@ -10,7 +10,8 @@
template <class Functional, class BitVector>
SolverFactory<Functional, BitVector>::SolverFactory(
const Dune::ParameterTree& parset,
- Functional& J) :
+ Functional& J,
+ const BitVector& ignoreNodes) :
J_(Dune::Solvers::wrap_own_share<const Functional>(std::forward<Functional>(J))),
//linear solver
//preconditioner_(),
@@ -24,6 +25,7 @@ SolverFactory<Functional, BitVector>::SolverFactory(
tnnmgStep_ = std::make_shared<Step>(*J_, dummyIterate_, nonlinearSmoother_, linearSolver_, DefectProjection(), LineSearchSolver());
tnnmgStep_->setPreSmoothingSteps(parset.get<int>("main.pre"));
+ tnnmgStep_->setIgnore(ignoreNodes);
}
template <class Functional, class BitVector>
diff --git a/src/spatial-solving/solverfactory.hh b/src/spatial-solving/solverfactory.hh
index c9a3d3b189c46d1a344b10419c7a05e3ac9b8008..ebf3034905b961fdface27db4e31b4cdeb774148 100644
--- a/src/spatial-solving/solverfactory.hh
+++ b/src/spatial-solving/solverfactory.hh
@@ -41,7 +41,8 @@ class SolverFactory {
SolverFactory(Dune::ParameterTree const &,
- Functional&);
+ Functional&,
+ const BitVector&);
void setProblem(Vector& x);
diff --git a/src/time-stepping/adaptivetimestepper.cc b/src/time-stepping/adaptivetimestepper.cc
index bfde1aebed7853c432d89d87c30a5146440d05bc..c4a7371e6ecf51f3a4f933ba6782d22be67b502d 100644
--- a/src/time-stepping/adaptivetimestepper.cc
+++ b/src/time-stepping/adaptivetimestepper.cc
@@ -21,6 +21,7 @@ template <class Factory, class NBodyAssembler, class Updaters, class ErrorNorm>
AdaptiveTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::AdaptiveTimeStepper(
Dune::ParameterTree const &parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
Updaters ¤t,
@@ -34,6 +35,7 @@ AdaptiveTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::AdaptiveTimeS
finalTime_(parset.get<double>("problem.finalTime")),
parset_(parset),
nBodyAssembler_(nBodyAssembler),
+ ignoreNodes_(ignoreNodes),
globalFriction_(globalFriction),
bodywiseNonmortarBoundaries_(bodywiseNonmortarBoundaries),
current_(current),
@@ -104,7 +106,7 @@ AdaptiveTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::step(
Updaters const &oldUpdaters, double rTime, double rTau) {
UpdatersWithCount newUpdatersAndCount = {oldUpdaters.clone(), {}};
newUpdatersAndCount.count =
- MyCoupledTimeStepper(finalTime_, parset_, nBodyAssembler_, globalFriction_, bodywiseNonmortarBoundaries_,
+ MyCoupledTimeStepper(finalTime_, parset_, nBodyAssembler_, ignoreNodes_, globalFriction_, bodywiseNonmortarBoundaries_,
newUpdatersAndCount.updaters, errorNorms_,
externalForces_)
.step(rTime, rTau);
diff --git a/src/time-stepping/adaptivetimestepper.hh b/src/time-stepping/adaptivetimestepper.hh
index d0a014d07c257e28f3f54fd4e52656e7d6f5bde3..f1f9b07ccb3b1dbb59b8deac5d182a42249d3989 100644
--- a/src/time-stepping/adaptivetimestepper.hh
+++ b/src/time-stepping/adaptivetimestepper.hh
@@ -23,6 +23,7 @@ class AdaptiveTimeStepper {
};
using Vector = typename Factory::Vector;
+ using IgnoreVector = typename Factory::BitVector;
//using ConvexProblem = typename Factory::ConvexProblem;
//using Nonlinearity = typename Factory::Nonlinearity;
@@ -36,6 +37,7 @@ class AdaptiveTimeStepper {
AdaptiveTimeStepper(
Dune::ParameterTree const &parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
Updaters ¤t,
@@ -58,6 +60,8 @@ class AdaptiveTimeStepper {
double finalTime_;
Dune::ParameterTree const &parset_;
NBodyAssembler& nBodyAssembler_;
+ const IgnoreVector& ignoreNodes_;
+
GlobalFriction& globalFriction_;
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries_;
diff --git a/src/time-stepping/coupledtimestepper.cc b/src/time-stepping/coupledtimestepper.cc
index 0ee30c60dd44ca16f32615321347458cfa174b3f..a102bb04f0c72c9a46239b7c5e3dffef101a5b6b 100644
--- a/src/time-stepping/coupledtimestepper.cc
+++ b/src/time-stepping/coupledtimestepper.cc
@@ -8,6 +8,7 @@ template <class Factory, class NBodyAssembler, class Updaters, class ErrorNorm>
CoupledTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::CoupledTimeStepper(
double finalTime, Dune::ParameterTree const &parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
Updaters updaters,
@@ -16,6 +17,7 @@ CoupledTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::CoupledTimeSte
: finalTime_(finalTime),
parset_(parset),
nBodyAssembler_(nBodyAssembler),
+ ignoreNodes_(ignoreNodes),
globalFriction_(globalFriction),
bodywiseNonmortarBoundaries_(bodywiseNonmortarBoundaries),
updaters_(updaters),
@@ -44,7 +46,7 @@ CoupledTimeStepper<Factory, NBodyAssembler, Updaters, ErrorNorm>::step(double re
updaters_.rate_->setup(ell, tau, newRelativeTime, velocityRHS, velocityIterate, velocityMatrix);
FixedPointIterator<Factory, NBodyAssembler, Updaters, ErrorNorm> fixedPointIterator(
- parset_, nBodyAssembler_, globalFriction_, bodywiseNonmortarBoundaries_, errorNorms_);
+ parset_, nBodyAssembler_, ignoreNodes_, globalFriction_, bodywiseNonmortarBoundaries_, errorNorms_);
auto const iterations = fixedPointIterator.run(updaters_, velocityMatrix,
velocityRHS, velocityIterate);
return iterations;
diff --git a/src/time-stepping/coupledtimestepper.hh b/src/time-stepping/coupledtimestepper.hh
index bc2b1f67c72655879b294bc5a1659c8d41c753c7..4610c7a77fe5c8558e5933fa9bb78a4bcab9b339 100644
--- a/src/time-stepping/coupledtimestepper.hh
+++ b/src/time-stepping/coupledtimestepper.hh
@@ -12,6 +12,7 @@ template <class Factory, class NBodyAssembler, class Updaters, class ErrorNorm>
class CoupledTimeStepper {
using Vector = typename Factory::Vector;
using Matrix = typename Factory::Matrix;
+ using IgnoreVector = typename Factory::BitVector;
//using Nonlinearity = typename Factory::Nonlinearity;
public:
using GlobalFriction = typename FixedPointIterator<Factory, NBodyAssembler, Updaters, ErrorNorm>::GlobalFriction;
@@ -22,6 +23,7 @@ class CoupledTimeStepper {
CoupledTimeStepper(double finalTime,
Dune::ParameterTree const &parset,
NBodyAssembler& nBodyAssembler,
+ const IgnoreVector& ignoreNodes,
GlobalFriction& globalFriction,
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries,
Updaters updaters, const std::vector<const ErrorNorm* >& errorNorms,
@@ -33,6 +35,8 @@ class CoupledTimeStepper {
double finalTime_;
Dune::ParameterTree const &parset_;
NBodyAssembler& nBodyAssembler_;
+ const IgnoreVector& ignoreNodes_;
+
GlobalFriction& globalFriction_;
const std::vector<const BitVector*>& bodywiseNonmortarBoundaries_;
diff --git a/src/time-stepping/state.cc b/src/time-stepping/state.cc
index 1e01e8a32519a9893591a86fd67565dc34cd612b..4184f3d765a8b014225ad146af0ca30683a7026e 100644
--- a/src/time-stepping/state.cc
+++ b/src/time-stepping/state.cc
@@ -26,9 +26,11 @@ auto initStateUpdater(
auto localUpdater = std::make_shared<AgeingLawStateUpdater<ScalarVector, Vector>>(
alpha_initial[nonmortarIdx], *contactCouplings[i]->nonmortarBoundary().getVertices(), coupling->frictionData().L(), coupling->frictionData().V0());
+ localUpdater->setBodyIndex(nonmortarIdx);
stateUpdater->addLocalUpdater(localUpdater);
}
+ break;
case Config::SlipLaw:
for (size_t i=0; i<couplings.size(); i++) {
const auto& coupling = couplings[i];
@@ -36,14 +38,17 @@ auto initStateUpdater(
auto localUpdater = std::make_shared<SlipLawStateUpdater<ScalarVector, Vector>>(
alpha_initial[nonmortarIdx], *contactCouplings[i]->nonmortarBoundary().getVertices(), coupling->frictionData().L(), coupling->frictionData().V0());
+ localUpdater->setBodyIndex(nonmortarIdx);
stateUpdater->addLocalUpdater(localUpdater);
}
+ break;
default:
assert(false);
+ break;
+ }
return stateUpdater;
- }
}
#include "state_tmpl.cc"
diff --git a/src/time-stepping/state/ageinglawstateupdater.cc b/src/time-stepping/state/ageinglawstateupdater.cc
index 6862029e930dac8ea8fe81fc6d918555a9b3fcac..a4c4737784142a11b7d4ce14393a9c728b9d740c 100644
--- a/src/time-stepping/state/ageinglawstateupdater.cc
+++ b/src/time-stepping/state/ageinglawstateupdater.cc
@@ -66,7 +66,9 @@ template <class ScalarVector, class Vector>
void AgeingLawStateUpdater<ScalarVector, Vector>::extractAlpha(
ScalarVector& alpha) {
- assert(alpha.size() == nodes_.size());
+ if (alpha.size() != nodes_.size()) {
+ alpha.resize(nodes_.size());
+ }
for (size_t i=0; i<localToGlobal_.size(); i++) {
alpha[localToGlobal_[i]] = alpha_[i];
diff --git a/src/time-stepping/state/stateupdater.hh b/src/time-stepping/state/stateupdater.hh
index ed7a5241b491df82e903f104a5ab494c63ff6196..458fc81a447aa833e27e10f7543439a1050ae887 100644
--- a/src/time-stepping/state/stateupdater.hh
+++ b/src/time-stepping/state/stateupdater.hh
@@ -57,7 +57,7 @@ template <class ScalarVectorTEMPLATE, class Vector> class StateUpdater {
}
}
- void extractAlpha(std::vector<ScalarVector>& alpha) {
+ void extractAlpha(std::vector<ScalarVector>& alpha) {
for (size_t i=0; i<localStateUpdaters_.size(); i++) {
auto& localStateUpdater = localStateUpdaters_[i];
localStateUpdater->extractAlpha(alpha[localStateUpdater->bodyIndex()]);