From b43335cad4d63c41328b80b2d2c8afecf5bc7c1a Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Tue, 25 Sep 2012 10:30:44 +0200
Subject: [PATCH] Remove Newmark step in the time stepping schem
This partially reverts d84ddb712891e2fe8a4bd067589a31c86e8a3eea
---
src/Makefile.am | 1 -
src/compute_state_dieterich_newmark.cc | 35 -------------------------
src/compute_state_dieterich_newmark.hh | 6 -----
src/compute_state_dieterich_newmark.tex | 34 ------------------------
src/one-body-sample.cc | 23 ++--------------
src/one-body-sample.parset | 1 -
6 files changed, 2 insertions(+), 98 deletions(-)
delete mode 100644 src/compute_state_dieterich_newmark.cc
delete mode 100644 src/compute_state_dieterich_newmark.hh
delete mode 100644 src/compute_state_dieterich_newmark.tex
diff --git a/src/Makefile.am b/src/Makefile.am
index ed547937..df019ae4 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -6,7 +6,6 @@ bin_PROGRAMS = \
SOURCES = \
assemblers.cc \
- compute_state_dieterich_newmark.cc \
compute_state_dieterich_euler.cc \
compute_state_dieterich_common.cc \
compute_state_ruina.cc \
diff --git a/src/compute_state_dieterich_newmark.cc b/src/compute_state_dieterich_newmark.cc
deleted file mode 100644
index 18704f48..00000000
--- a/src/compute_state_dieterich_newmark.cc
+++ /dev/null
@@ -1,35 +0,0 @@
-#include <cassert>
-
-#include <dune/tnnmg/problem-classes/bisection.hh>
-
-#include <dune/tectonic/mydirectionalconvexfunction.hh>
-
-#include "compute_state_dieterich_common.hh"
-#include "compute_state_dieterich_newmark.hh"
-
-double state_update_dieterich_newmark_bisection(double tau, double VoL,
- double old_state,
- double old_stated) {
- DieterichNonlinearity::VectorType const start(0);
- DieterichNonlinearity::VectorType const direction(1);
- DieterichNonlinearity const phi(VoL);
- MyDirectionalConvexFunction<DieterichNonlinearity> const J(
- 2.0 / tau, old_stated + 2.0 / tau * old_state, phi, start, direction);
-
- int bisectionsteps = 0;
- Bisection const bisection(0.0, 1.0, 1e-12, true, 0); // TODO
- return bisection.minimize(J, 0.0, 0.0, bisectionsteps); // TODO
-}
-
-double state_update_dieterich_newmark(double tau, double VoL, double old_state,
- double old_stated) {
- double const new_state =
- state_update_dieterich_newmark_bisection(tau, VoL, old_state, old_stated);
- double const new_stated = 2.0 / tau * (new_state - old_state) - old_stated;
- double const rhs = -new_stated;
- double const lhs = VoL - std::exp(-new_state);
-
- assert(std::abs(lhs - rhs) < 1e-10);
-
- return new_state;
-}
diff --git a/src/compute_state_dieterich_newmark.hh b/src/compute_state_dieterich_newmark.hh
deleted file mode 100644
index 755dbf03..00000000
--- a/src/compute_state_dieterich_newmark.hh
+++ /dev/null
@@ -1,6 +0,0 @@
-#ifndef COMPUTE_STATE_DIETERICH_NEWMARK_HH
-#define COMPUTE_STATE_DIETERICH_NEWMARK_HH
-
-double state_update_dieterich_newmark(double h, double VoL, double old_state,
- double old_stated);
-#endif
diff --git a/src/compute_state_dieterich_newmark.tex b/src/compute_state_dieterich_newmark.tex
deleted file mode 100644
index 242721bc..00000000
--- a/src/compute_state_dieterich_newmark.tex
+++ /dev/null
@@ -1,34 +0,0 @@
-\documentclass{scrartcl}
-\usepackage{amsmath}
-\begin{document}
-\noindent The Newmark scheme in its classical form with $\gamma = 1/2$
-and $\beta = 1/4$ reads
-\begin{align}
- \label{eq:1} \dot \alpha_1
- &= \dot \alpha_0 + \frac \tau 2 (\ddot \alpha_0 + \ddot \alpha_1 )\\
- \label{eq:2} \alpha_1
- &= \alpha_0 + \tau \dot \alpha_0 + \frac {\tau^2}4 ( \ddot \alpha_0 + \ddot \alpha_1 )\text.
- \intertext{We can also write \eqref{eq:2} as}
- \nonumber \ddot \alpha_1 + \ddot \alpha_0
- &= \frac 4{\tau^2} ( \alpha_1 - \alpha_0 - \tau \dot \alpha_0)
- \intertext{so that it yields}
- \label{eq:3} \dot \alpha_1
- &= \dot \alpha_0 + \frac 2\tau ( \alpha_1 - \alpha_0 - \tau \dot \alpha_0) = \frac 2\tau ( \alpha_1 - \alpha_0) - \dot \alpha_0
-\end{align}
-in conjunction with \eqref{eq:1}. The problem
-\begin{align*}
- -\dot \alpha_1 \in \partial \psi(\alpha_1)
-\end{align*}
-then becomes
-\begin{align*}
- \dot \alpha_0 -\frac 2\tau ( \alpha_1 - \alpha_0)
- &\in \partial \psi(\alpha_1)\\
- \psi(\beta) - \psi(\alpha_1)
- &\ge (\dot \alpha_0 -\frac 2\tau ( \alpha_1 - \alpha_0), \beta - \alpha_1)
- \quad \forall \beta\\
- \frac 2\tau ( \alpha_1, \beta - \alpha_1) + \psi(\beta) - \psi(\alpha_1)
- &\ge (\dot \alpha_0 + \frac 2\tau \alpha_0, \beta - \alpha_1)
- \quad \forall \beta
-\end{align*}
-After which $\dot \alpha_1$ can be computed according to \eqref{eq:3}.
-\end{document}
diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index 8551d62e..b3d1420e 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -66,7 +66,6 @@
#include "assemblers.hh"
#include "compute_state_dieterich_euler.hh"
-#include "compute_state_dieterich_newmark.hh"
#include "compute_state_ruina.hh"
#include "mysolver.hh"
#include "vtk.hh"
@@ -256,10 +255,6 @@ int main(int argc, char *argv[]) {
alpha_old = parset.get<double>("boundary.friction.state.initial");
SingletonVectorType alpha(alpha_old);
- SingletonVectorType alphad_old(finestSize);
- alphad_old = 0.0;
- SingletonVectorType alphad(alphad_old);
-
SingletonVectorType vonMisesStress;
// }}}
@@ -380,21 +375,8 @@ int main(int argc, char *argv[]) {
switch (parset.get<Config::state_model>(
"boundary.friction.state.model")) {
case Config::Dieterich:
- switch (
- parset.get<Config::scheme>("stateTimeSteppingScheme")) {
- case Config::ImplicitEuler:
- alpha[i] = state_update_dieterich_euler(tau, V / L,
- alpha_old[i]);
- break;
- case Config::Newmark:
- alpha[i] = state_update_dieterich_newmark(
- tau, V / L, alpha_old[i], alphad_old[i]);
- alphad[i] =
- 2.0 / tau * (alpha[i] - alpha_old[i]) - alphad_old[i];
- break;
- default:
- assert(false); // Nothing else implemented
- }
+ alpha[i] =
+ state_update_dieterich_euler(tau, V / L, alpha_old[i]);
break;
case Config::Ruina:
alpha[i] = state_update_ruina(tau, V * tau / L, alpha_old[i]);
@@ -449,7 +431,6 @@ int main(int argc, char *argv[]) {
}
alpha_old = alpha;
- alphad_old = alphad;
{ // Write all kinds of data
for (size_t i = 0; i < frictionalNodes.size(); ++i)
diff --git a/src/one-body-sample.parset b/src/one-body-sample.parset
index 09740d21..8109f404 100644
--- a/src/one-body-sample.parset
+++ b/src/one-body-sample.parset
@@ -11,7 +11,6 @@ printVelocitySteppingComparison = false
enableTimer = false
timeSteppingScheme = implicitEuler
-stateTimeSteppingScheme = implicitEuler
[grid]
refinements = 4
--
GitLab