A bit of text on the dynamic model

doc/relaxed-micromorphic-continuum.bib

@@ -8,3 +8,11 @@
  pages = {53--84},
  doi = {10.1093/qjmam/hbu027}
 }
+
+@Book{hairer_lubich_wanner:2006,
+ title = {Geometric Numerical Integration---Structure-Preserving Algorithms for Ordinary Differential Equations},
+ publisher = {Springer},
+ year = {2016},
+ author = {Ernst Hairer and Christian Lubich and Gerhard Wanner},
+ edition = {Zweite Auflage},
+}

doc/relaxed-micromorphic-continuum.tex

@@ -88,10 +88,9 @@ matrix curl as
  (\Curl P)_i
  \colonequals
  \Big( \frac{\partial P_{i3}}{\partial x_2} - \frac{\partial P_{i2}}{\partial x_3}, \quad
-       \frac{\partial P_{i3}}{\partial x_2} - \frac{\partial P_{i3}}{\partial x_2}, \quad
-       \frac{\partial P_{i3}}{\partial x_2} - \frac{\partial P_{i3}}{\partial x_2} \Big).
+       \frac{\partial P_{i1}}{\partial x_3} - \frac{\partial P_{i3}}{\partial x_1}, \quad
+       \frac{\partial P_{i2}}{\partial x_1} - \frac{\partial P_{i1}}{\partial x_2} \Big).
 \end{equation*}
-\todosander{Wie würde das Modell in 2d aussehen?}
 
 We define the quantities
 \begin{alignat*}{2}
@@ -224,5 +223,138 @@ norm $\tnorm{\cdot}_\mathcal{X}$.
 
 \section{The dynamic case}
 
+We now investigate a time-dependent relaxed micromorphic continuum, using the
+Hamilton formalism.  As shown above, the model is described by a quadratic potential energy,
+and we assume the kinetic energy to be quadratic as well.
+
+\subsection{The implicit midpoint rule for quadratic Hamiltonians}
+
+Consider a model with quadratic potential energy
+\begin{equation*}
+ U(q) = \frac{1}{2} q^TA q + b^T q
+\end{equation*}
+and quadratic kinetic energy
+\begin{equation*}
+ T(q, \dot{q}) = T(\dot{q}) = \frac{1}{2} \dot{q}^T M \dot{q},
+\end{equation*}
+where $A$ and $M$ are symmetric and positive definite matrices.
+
+As the kinetic energy is quadratic, the corresponding Hamiltonian
+is simply the sum of the two energies~\cite[Chapter~VI.1.2]{hairer_lubich_wanner:2006}
+Using that $p = M \dot{q}$ we get
+\begin{equation*}
+ H(p,q)
+ =
+ T(\dot{q}(p)) + U(q)
+ =
+ \frac{1}{2} p^TM^{-T} p + \frac{1}{2} q^T A q + b^T q.
+\end{equation*}
+
+We want to integrate this system using an implicit midpoint rule.  As is
+well known, such an integrator is symplectic, second-order accurate~\cite[Theorem~VI.3.5]{hairer_lubich_wanner:2006}
+and B-stable~\cite{}.  One step of this method is
+\begin{equation*}
+ \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix}
+ =
+ \begin{pmatrix} p^k \\ q^k \end{pmatrix}
+ +
+ \tau J^{-1} \nabla H\Big(\frac{p^k + p^{k+1}}{2}, \frac{q^k + q^{k+1}}{2}\Big),
+\end{equation*}
+where $\tau$ is the time step size, and $J^{-1}$ is the matrix
+\begin{equation*}
+ J^{-1}
+ \colonequals
+ \begin{pmatrix} 0 & -I \\ I & 0 \end{pmatrix}.
+\end{equation*}
+With
+\begin{equation*}
+ \nabla H(p,q)
+ =
+ \begin{pmatrix}
+  M^{-T} p \\ Aq + b
+ \end{pmatrix}
+\end{equation*}
+the integrator becomes
+\begin{equation*}
+ \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix}
+ - \frac{\tau}{2} \begin{pmatrix} 0 & -I \\ I & 0 \end{pmatrix}
+ \underbrace{
+ \begin{pmatrix}
+  M^{-T} p^{k+1} \\ Aq^{k+1}
+ \end{pmatrix}
+ }_{ = \begin{pmatrix} M^{-T} & 0 \\ 0 & A \end{pmatrix} \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix} }
+ =
+ \begin{pmatrix} p^k \\ q^k \end{pmatrix}
+ + \frac{\tau}{2} \begin{pmatrix} 0 & -I \\ I & 0 \end{pmatrix}
+ \begin{pmatrix}
+  M^{-T} p^k \\ Aq^k + 2b
+ \end{pmatrix}
+\end{equation*}
+This can be transformed to
+\begin{equation*}
+ \Bigg[
+ \begin{pmatrix} I & 0 \\ 0 & I \end{pmatrix}
+ - \frac{\tau}{2}
+ \begin{pmatrix}
+  0 & -A \\ M^{-T} & 0
+ \end{pmatrix}
+ \Bigg]
+ \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix}
+ =
+ \begin{pmatrix}
+  p^k  - \frac{\tau}{2} A q^k - \tau b\\
+  q^k + \frac{\tau}{2} M^{-T} p^k
+ \end{pmatrix},
+\end{equation*}
+and this in turn to
+\begin{equation*}
+ \begin{pmatrix}
+   I & \frac{\tau}{2} A \\
+   - \frac{\tau}{2} M^{-T} & I
+ \end{pmatrix}
+ \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix}
+ =
+ \begin{pmatrix}
+  p^k  - \frac{\tau}{2} A q^k - \tau b\\
+  q^k + \frac{\tau}{2} M^{-T} p^k
+ \end{pmatrix},
+\end{equation*}
+Eliminate the lower left block:
+\begin{equation*}
+ \begin{pmatrix}
+   I & \frac{\tau}{2} A \\
+   0 & I + \frac{\tau^2}{4} M^{-T}A
+ \end{pmatrix}
+ \begin{pmatrix} p^{k+1} \\ q^{k+1} \end{pmatrix}
+ =
+ \begin{pmatrix}
+  p^k  - \frac{\tau}{2} A q^k - \tau b\\
+  q^k + \frac{\tau}{2} M^{-T} p^k + \frac{\tau}{2} M^{-T}(p^k - \frac{\tau}{2}Aq^k - \tau b)
+ \end{pmatrix}.
+\end{equation*}
+
+This can be solved by backward substitution.  In algorithmic form this reads:
+\begin{enumerate}
+ \item Solve
+  \begin{equation*}
+   \Big(M^T + \frac{\tau^2}{4} A \Big) q^{k+1}
+   =
+   M^Tq^k + \tau p^k - \frac{\tau^2}{4}Aq^k - \frac{\tau^2}{2}b
+  \end{equation*}
+
+ \item Set
+  \begin{equation*}
+   p^{k+1}
+   =
+   p^k - \frac{\tau}{2} Aq^k - 2b - \frac{\tau}{2}Aq^{k+1}.
+  \end{equation*}
+
+\end{enumerate}
+
+
+
+
+
+
 \printbibliography
 \end{document}