diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index a815fb4b38573cc977d005ba4e83595796a216d9..9fc19e0b9b338e355a36f964adae56be3dd20585 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -189,20 +189,20 @@ int main(int argc, char *argv[]) {
leafView, p1Basis, frictionalNodes);
// {{{ Initialise vectors
- VectorType u4(finestSize);
- u4 = 0.0; // Has to be zero!
+ VectorType u(finestSize);
+ u = 0.0; // Has to be zero!
- SingletonVectorType s4_old(finestSize);
- s4_old = parset.get<double>("boundary.friction.state.initial");
+ SingletonVectorType alpha_old(finestSize);
+ alpha_old = parset.get<double>("boundary.friction.state.initial");
- VectorType u4_diff(finestSize);
- u4_diff = 0.0; // Has to be zero!
+ VectorType u_diff(finestSize);
+ u_diff = 0.0; // Has to be zero!
- auto s4_new = Dune::make_shared<SingletonVectorType>(finestSize);
- *s4_new = s4_old;
+ auto alpha_new = Dune::make_shared<SingletonVectorType>(finestSize);
+ *alpha_new = alpha_old;
SingletonVectorType vonMisesStress;
- VectorType b4;
+ VectorType rhs;
// }}}
typedef MyConvexProblem<MatrixType, VectorType> MyConvexProblemType;
@@ -253,13 +253,13 @@ int main(int argc, char *argv[]) {
{
assemble_neumann<GridType, GridView, SmallVector, P1Basis>(
- leafView, p1Basis, neumannNodes, b4, neumannFunction, time);
- stiffnessMatrix.mmv(u4, b4);
+ leafView, p1Basis, neumannNodes, rhs, neumannFunction, time);
+ stiffnessMatrix.mmv(u, rhs);
// Apply Dirichlet condition
for (size_t i = 0; i < finestSize; ++i)
if (ignoreNodes[i].count() == dim) {
- dirichletFunction.evaluate(time, u4_diff[i][0]);
- u4_diff[i][0] /= refinement_factor;
+ dirichletFunction.evaluate(time, u_diff[i][0]);
+ u_diff[i][0] /= refinement_factor;
}
for (size_t state_fpi = 0;
@@ -268,15 +268,16 @@ int main(int argc, char *argv[]) {
++state_fpi) {
auto myGlobalNonlinearity =
assemble_nonlinearity<MatrixType, VectorType>(
- parset.sub("boundary.friction"), nodalIntegrals, s4_new, h);
+ parset.sub("boundary.friction"), nodalIntegrals, alpha_new,
+ h);
MyConvexProblemType const myConvexProblem(stiffnessMatrix,
- *myGlobalNonlinearity, b4);
+ *myGlobalNonlinearity, rhs);
MyBlockProblemType myBlockProblem(parset, myConvexProblem);
auto multigridStep = mySolver.getSolver();
- multigridStep->setProblem(u4_diff, myBlockProblem);
+ multigridStep->setProblem(u_diff, myBlockProblem);
- VectorType const u4_diff_old = u4_diff;
+ VectorType const u_diff_old = u_diff;
LoopSolver<VectorType> overallSolver(
multigridStep, parset.get<size_t>("solver.tnnmg.maxiterations"),
solver_tolerance, &energyNorm, verbosity,
@@ -285,7 +286,7 @@ int main(int argc, char *argv[]) {
for (size_t i = 0; i < frictionalNodes.size(); ++i) {
if (frictionalNodes[i][0]) {
- double const unorm = u4_diff[i].two_norm();
+ double const unorm = u_diff[i].two_norm();
// // the (logarithmic) steady state corresponding to the
// // current velocity
@@ -294,11 +295,12 @@ int main(int argc, char *argv[]) {
switch (parset.get<Config::state_model>(
"boundary.friction.state.model")) {
case Config::Dieterich:
- (*s4_new)[i] =
- state_update_dieterich(h, unorm / L, s4_old[i]);
+ (*alpha_new)[i] =
+ state_update_dieterich(h, unorm / L, alpha_old[i]);
break;
case Config::Ruina:
- (*s4_new)[i] = state_update_ruina(h, unorm / L, s4_old[i]);
+ (*alpha_new)[i] =
+ state_update_ruina(h, unorm / L, alpha_old[i]);
break;
}
}
@@ -307,7 +309,7 @@ int main(int argc, char *argv[]) {
std::cout << '.';
std::cout.flush();
}
- if (energyNorm.diff(u4_diff_old, u4_diff) <
+ if (energyNorm.diff(u_diff_old, u_diff) <
parset.get<double>("solver.tnnmg.fixed_point_tolerance"))
break;
}
@@ -317,8 +319,8 @@ int main(int argc, char *argv[]) {
if (parset.get<bool>("writeEvolution")) {
double out;
neumannFunction.evaluate(time, out);
- octave_writer << (*s4_new)[first_frictional_node][0] << " "
- << u4[first_frictional_node][0] * 1e6 << " " << out
+ octave_writer << (*alpha_new)[first_frictional_node][0] << " "
+ << u[first_frictional_node][0] * 1e6 << " " << out
<< std::endl;
}
@@ -326,9 +328,9 @@ int main(int argc, char *argv[]) {
// with the Ruina state evolution law.
// Jumps at 120 and 360 timesteps; v1 = .6 * v2;
if (parset.get<bool>("printVelocitySteppingComparison")) {
- double const v = u4_diff[first_frictional_node].two_norm() / h / L;
+ double const v = u_diff[first_frictional_node].two_norm() / h / L;
- double const euler = (*s4_new)[first_frictional_node];
+ double const euler = (*alpha_new)[first_frictional_node];
double direct;
if (run < 120) {
direct = euler;
@@ -350,28 +352,28 @@ int main(int argc, char *argv[]) {
// Record the coefficient of friction at a fixed node
if (parset.get<bool>("printCoefficient")) {
- double const V = u4_diff[first_frictional_node].two_norm();
- double const state = (*s4_new)[first_frictional_node];
+ double const V = u_diff[first_frictional_node].two_norm();
+ double const state = (*alpha_new)[first_frictional_node];
coefficient_writer << (mu + a * std::log(V * eta) +
b * (state - std::log(eta * L))) << std::endl;
}
}
- u4 += u4_diff;
- s4_old = *s4_new;
+ u += u_diff;
+ alpha_old = *alpha_new;
// Compute von Mises stress and write everything to a file
if (parset.get<bool>("writeVTK")) {
VonMisesStressAssembler<GridType> localStressAssembler(
E, nu,
Dune::make_shared<BasisGridFunction<P1Basis, VectorType> const>(
- p1Basis, u4));
+ p1Basis, u));
FunctionalAssembler<P0Basis>(p0Basis)
.assemble(localStressAssembler, vonMisesStress, true);
writeVtk<P1Basis, P0Basis, VectorType, SingletonVectorType, GridView>(
- p1Basis, u4, *s4_new, p0Basis, vonMisesStress, leafView,
+ p1Basis, u, *alpha_new, p0Basis, vonMisesStress, leafView,
(boost::format("obs%d") % run).str());
}
}
@@ -381,10 +383,10 @@ int main(int argc, char *argv[]) {
if (parset.get<bool>("printFrictionalBoundary")) {
// Print displacement on frictional boundary
- boost::format const formatter("u4[%02d] = %+3e");
+ boost::format const formatter("u[%02d] = %+3e");
for (size_t i = 0; i < frictionalNodes.size(); ++i)
if (frictionalNodes[i][0])
- std::cout << (boost::format(formatter) % i % u4[i]) << std::endl;
+ std::cout << (boost::format(formatter) % i % u[i]) << std::endl;
}
octave_writer.close();