diff --git a/doc/manual/manual-appendix-materials-cohesive.tex b/doc/manual/manual-appendix-materials-cohesive.tex
new file mode 100644
index 000000000..1cf1f3aae
--- /dev/null
+++ b/doc/manual/manual-appendix-materials-cohesive.tex
@@ -0,0 +1,30 @@
+\section{Cohesive linear}
+
+\begin{MaterialDesc}{cohesive\_linear}{ssect:smm:cl:coh-snozzi}
+\matparam{sigma\_c}{Real}{Critical stress}
+\matparam{beta}{Real}{$\beta$ parameter}
+\matparam{G\_cI}{Real}{Mode I fracture energy}
+\matparam{G\_cII}{Real}{Mode II fracture energy}
+\matparam{kappa}{Real}{$\kappa$ parameter}
+\matparam{penalty}{Real}{penalty coefficient $\alpha$}
+\end{MaterialDesc}
+
+\section{Cohesive bilinear}
+
+\begin{MaterialDesc}{cohesive\_bilinear}{ssect:smm:cl:coh-snozzi}
+\matparam{sigma\_c}{Real}{Critical stress}
+\matparam{beta}{Real}{$\beta$ parameter}
+\matparam{G\_cI}{Real}{Mode I fracture energy}
+\matparam{G\_cII}{Real}{Mode II fracture energy}
+\matparam{kappa}{Real}{$\kappa$ parameter}
+\matparam{penalty}{Real}{Penalty coefficient $\alpha$}
+\matparam{delta\_0}{Real}{Elastic limit displacement $\delta_0$}
+\end{MaterialDesc}
+
+\section{Cohesive exponential}
+
+\begin{MaterialDesc}{cohesive\_exponential}{ssect:smm:cl:coh-exponential}
+\matparam{sigma\_c}{Real}{Critical stress}
+\matparam{beta}{Real}{$\beta$ parameter}
+\matparam{delta\_c}{Real}{Critical displacement $\delta_\mathrm{c}$}
+\end{MaterialDesc}
diff --git a/doc/manual/manual-appendix-materials-extra-materials.tex b/doc/manual/manual-appendix-materials-extra-materials.tex
new file mode 100644
index 000000000..bf93b615d
--- /dev/null
+++ b/doc/manual/manual-appendix-materials-extra-materials.tex
@@ -0,0 +1,41 @@
+\section{Linear elastic orthotropic}
+
+\begin{MaterialDesc}{elastic\_orthotropic}{}
+\matparam{rho}{Real}{Density}
+\matparam{n1}{Real}{Direction of main material axis}
+\matparam{n2}{Real}{Direction of secondary material axis}
+\matparam{n3}{Real}{Direction of tertiery material axis}
+\matparam{alpha}{Real}{Proportion of viscous stress}
+\matparam{E1}{Real}{Young's modulus $E_1$}
+\matparam{E2}{Real}{Young's modulus $E_2$}
+\matparam{E3}{Real}{Young's modulus $E_3$}
+\matparam{nu12}{Real}{Poisson's ration $\nu_{12}$}
+\matparam{nu13}{Real}{Poisson's ration $\nu_{13}$}
+\matparam{nu23}{Real}{Poisson's ration $\nu_{23}$}
+\matparam{G12}{Real}{Shear modulus $G_{12}$}
+\matparam{G13}{Real}{Shear modulus $G_{13}$}
+\matparam{G23}{Real}{Shear modulus $G_{23}$}
+\end{MaterialDesc}
+
+\section{Linear elastic anisotropic}
+
+\begin{MaterialDesc}{elastic\_anisotropic}{}
+\matparam{rho}{Real}{Density}
+\matparam{n1}{Real}{Direction of main material axis}
+\matparam{n2}{Real}{Direction of secondary material axis}
+\matparam{n3}{Real}{Direction of tertiery material axis}
+\matparam{alpha}{Real}{Proportion of viscous stress}
+\matparam{Cij}{Real}{Stiffness matrix coefficients $C_{ij}$ ($i,j = 1,2,$ \dots voigt size)}
+\end{MaterialDesc}
+
+\section{Visco-plastic}
+
+\begin{MaterialDesc}{visco\_plastic}{}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{Real}{Plane stress simplification (only 2D problems)}
+\matparam{rate}{Real}{Rate sensitivity component}
+\matparam{edot0}{Real}{Reference strain rate}
+\matparam{ts}{Real}{Time}
+\end{MaterialDesc}
diff --git a/doc/manual/manual-appendix-materials-non-local.tex b/doc/manual/manual-appendix-materials-non-local.tex
new file mode 100644
index 000000000..2da1b8d24
--- /dev/null
+++ b/doc/manual/manual-appendix-materials-non-local.tex
@@ -0,0 +1,11 @@
+\section{Non-Local damage: Marigo}
+
+\begin{MaterialDesc}{marigo_non_local}{ssect:smm:cl:damage-marigo-non-local}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\matparam{Yd}{Random}{Rupture criterion}
+\matparam{S}{Real}{Damage Energy}
+\matparam{weight_function}{Subsection}{specified by the optional flavor}
+\end{MaterialDesc}
diff --git a/doc/manual/manual-appendix-materials.tex b/doc/manual/manual-appendix-materials.tex
index 4b7679bec..f27c11633 100644
--- a/doc/manual/manual-appendix-materials.tex
+++ b/doc/manual/manual-appendix-materials.tex
@@ -1,145 +1,80 @@
\chapter{Material parameters}
\label{app:material-parameters}
\section{Linear elastic isotropic}
-Keyword: \code{elastic}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{E}: Young's modulus
-\item \code{nu}: Poisson's ratio
-\item \code{Plane\_Stress}: Plane stress simplification (only 2D problems)
-\end{itemize}
-
-\section{Linear elastic orthotropic}
-
-Keyword: \code{elastic\_orthotropic}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{n1}: Direction of main material axis
-\item \code{n2}: Direction of secondary material axis
-\item \code{n3}: Direction of tertiery material axis
-\item \code{alpha}: Proportion of viscous stress
-\item \code{E1}: Young's modulus $E_1$
-\item \code{E2}: Young's modulus $E_2$
-\item \code{E3}: Young's modulus $E_3$
-\item \code{nu12}: Poisson's ration $\nu_{12}$
-\item \code{nu13}: Poisson's ration $\nu_{13}$
-\item \code{nu23}: Poisson's ration $\nu_{23}$
-\item \code{G12}: Shear modulus $G_{12}$
-\item \code{G13}: Shear modulus $G_{13}$
-\item \code{G23}: Shear modulus $G_{23}$
-\end{itemize}
-
-\section{Linear elastic anisotropic}
-
-Keyword: \code{elastic\_anisotropic}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{n1}: Direction of main material axis
-\item \code{n2}: Direction of secondary material axis
-\item \code{n3}: Direction of tertiery material axis
-\item \code{alpha}: Proportion of viscous stress
-\item \code{Cij}: Stiffness matrix coefficients $C_{ij}$ ($i,j = 1,2,$ \dots voigt size)
-\end{itemize}
+\begin{MaterialDesc}{elastic}{ssect:smm:linear-elastic}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\end{MaterialDesc}
-\section{Neohookean (finite strains)}
-Keyword: \code{neohookean}\\
+\section{Neohookean (finite strains)}
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{E}: Young's modulus
-\item \code{nu}: Poisson's ratio
-\item \code{Plane\_Stress}: Plane stress simplification (only 2D problems)
-\end{itemize}
+\begin{MaterialDesc}{neohookean}{ssect:smm:cl:neohookean}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\end{MaterialDesc}
\section{Standard linear solid}
-Keyword: \code{sls\_deviatoric}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{E}: Young's modulus
-\item \code{nu}: Poisson's ratio
-\item \code{Plane\_Stress}: Plane stress simplification (only 2D problems)
-\item \code{Eta}: Viscosity
-\item \code{Ev}: Stiffness of the viscous element
-\end{itemize}
+\begin{MaterialDesc}{sls\_deviatoric}{ssect:smm:cl:sls}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\matparam{Eta}{Real}{Viscosity}
+\matparam{Ev}{Real}{Stiffness of the viscous element}
+\end{MaterialDesc}
\section{Elasto-plastic linear isotropic hardening}
-Keyword: \code{plastic\_linear\_isotropic\_hardening}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{E}: Young's modulus
-\item \code{nu}: Poisson's ratio
-\item \code{Plane\_Stress}: Plane stress simplification (only 2D problems)
-\item \code{h}: Hardening modulus
-\item \code{sigma\_y}: Yielding stress
-\end{itemize}
-
-\section{Visco-plastic}
-
-Keyword: \code{visco\_plastic}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{rho}: Density
-\item \code{E}: Young's modulus
-\item \code{nu}: Poisson's ratio
-\item \code{Plane\_Stress}: Plane stress simplification (only 2D problems)
-\item \code{rate}: Rate sensitivity component
-\item \code{edot0}: Reference strain rate
-\item \code{ts}: Time
-\end{itemize}
-
-\section{Cohesive linear}
-
-Keyword: \code{cohesive\_linear}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{sigma\_c}: Critical stress
-\item \code{beta}: $\beta$ parameter
-\item \code{G\_cI}: Mode I fracture energy
-\item \code{G\_cII}: Mode II fracture energy
-\item \code{kappa}: $\kappa$ parameter
-\item \code{penalty}: penalty coefficient $\alpha$
-\end{itemize}
-
-\section{Cohesive bilinear}
-
-Keyword: \code{cohesive\_bilinear}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{sigma\_c}: Critical stress
-\item \code{beta}: $\beta$ parameter
-\item \code{G\_cI}: Mode I fracture energy
-\item \code{G\_cII}: Mode II fracture energy
-\item \code{kappa}: $\kappa$ parameter
-\item \code{penalty}: Penalty coefficient $\alpha$
-\item \code{delta\_0}: Elastic limit displacement $\delta_0$
-\end{itemize}
-
-\section{Cohesive exponential}
-
-Keyword: \code{cohesive\_exponential}\\
-
-\noindent Parameters:
-\begin{itemize}
-\item \code{sigma\_c}: Critical stress
-\item \code{beta}: $\beta$ parameter
-\item \code{delta\_c}: Critical displacement $\delta_\mathrm{c}$
-\end{itemize}
+\begin{MaterialDesc}{plastic\_linear\_isotropic\_hardening}{ssect:smm:cl:plastic}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\matparam{h}{Real}{Hardening modulus}
+\matparam{sigma\_y}{Real}{Yielding stress}
+\end{MaterialDesc}
+
+\section{Damage: Marigo}
+
+\begin{MaterialDesc}{marigo}{ssect:smm:cl:damage-marigo}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{Plane\_Stress}{bool}{Plane stress simplification (only 2D problems)}
+\matparam{Yd}{Random}{Rupture criterion}
+\matparam{S}{Real}{Damage Energy}
+\end{MaterialDesc}
+
+
+
+\section{Damage: Mazars}
+
+\begin{MaterialDesc}{mazars}{ssect:smm:cl:damage-mazars}
+\matparam{rho}{Real}{Density}
+\matparam{E}{Real}{Young's modulus}
+\matparam{nu}{Real}{Poisson's ratio}
+\matparam{At}{Real}{Traction post-peak asymptotic value}
+\matparam{Bt}{Real}{Traction decay shape}
+\matparam{Ac}{Real}{Compression post-peak asymptotic value}
+\matparam{Bc}{Real}{Compression decay shape}
+\matparam{K0}{Real}{Damage threshold}
+\matparam{beta}{Real}{Shear parameter}
+\end{MaterialDesc}
+
+\IfFileExists{manual-appendix-materials-extra-materials.tex}{\input{manual-appendix-materials-extra-materials}}{}
+
+\IfFileExists{manual-appendix-materials-cohesive.tex}{\input{manual-appendix-materials-cohesive}}{}
+
+
+%%% Local Variables:
+%%% mode: latex
+%%% TeX-master: "manual"
+%%% End:
diff --git a/doc/manual/manual-cohesive_laws.tex b/doc/manual/manual-cohesive_laws.tex
index 8a778d771..5049de84e 100644
--- a/doc/manual/manual-cohesive_laws.tex
+++ b/doc/manual/manual-cohesive_laws.tex
@@ -1,140 +1,140 @@
\subsection{Cohesive laws}
\label{sec:cohesive-laws}
-\subsubsection{Snozzi-Molinari Law}
+\subsubsection{Snozzi-Molinari Law}\label{ssect:smm:cl:coh-snozzi}
\begin{figure}[!hbt]
\centering
\subfloat[Linear]{\includegraphics[width=0.4\textwidth]{figures/linear_cohesive_law}}
\qquad
\subfloat[Bilinear]{\includegraphics[width=0.4\textwidth]{figures/bilinear_cohesive_law}}
\caption{Irreversible cohesive laws for explicit simulations.}
\label{fig:smm:coh:linear_cohesive_law}
\end{figure}
\akantu includes the Snozzi-Molinari~\cite{snozzi_cohesive_2013}
linear irreversible cohesive law (see
Figure~\ref{fig:smm:coh:linear_cohesive_law}). It is an extension to
the Camacho-Ortiz~\cite{camacho_computational_1996} cohesive law in
order to make dissipated fracture energy path-dependent. The concept
of free potential energy is dropped and a new independent parameter
$\kappa$ is introduced:
\begin{equation}
\kappa = \frac{G_\mathrm{c, II}}{G_\mathrm{c, I}}
\end{equation}
where $G_\mathrm{c, I}$ and $G_\mathrm{c, II}$ are respectively the
necessary works of separation per unit area to open completely a
cohesive zone under mode I and mode II. Their model yields to the
following equation for cohesive tractions $\vec{T}$ in case of crack
opening ${\delta}$:
\begin{equation}
\label{eq:smm:coh:tractions}
\vec{T} = \left( \frac{\beta^2}{\kappa} \Delta_\mathrm{t} \vec{t} +
\Delta_\mathrm{n} \vec{n} \right)
\frac{\sigma_\mathrm{c}}{\delta}
\left( 1- \frac{\delta}{\delta_\mathrm{c}} \right)
= \hat{\vec{T}}\, \frac{\sigma_\mathrm{c}}{\delta}
\left( 1- \frac{\delta}{\delta_\mathrm{c}} \right)
\end{equation}
where $\sigma_\mathrm{c}$ is the material strength along the fracture,
$\delta_\mathrm{c}$ the critical effective displacement after which
cohesive tractions are zero (complete decohesion), $\Delta_\mathrm{t}$
and $\Delta_\mathrm{n}$ are respectively the tangential and normal
components of the opening displacement vector $\vec{\Delta}$. $\beta$
is a weight that indicates how big the tangential opening contribution
is. The effective opening displacement is:
\begin{equation}
\delta = \sqrt{\frac{\beta^2}{\kappa^2} \Delta_\mathrm{t}^2 +
\Delta_\mathrm{n}^2}
\end{equation}
In case of unloading or reloading $\delta < \delta_\mathrm{max}$,
tractions are calculated as:
\begin{align}
T_\mathrm{n} &= \Delta_\mathrm{n}\,
\frac{\sigma_\mathrm{c}}{\delta_\mathrm{max}}
\left( 1- \frac{\delta_\mathrm{max}}{\delta_\mathrm{c}} \right) \\
T_\mathrm{t} &= \frac{\beta^2}{\kappa}\, \Delta_\mathrm{t}\,
\frac{\sigma_\mathrm{c}}{\delta_\mathrm{max}}
\left( 1- \frac{\delta_\mathrm{max}}{\delta_\mathrm{c}} \right)
\end{align}
so that they vary linearly between the origin and the maximum attained
tractions. As shown in Figure~\ref{fig:smm:coh:linear_cohesive_law},
in this law, the dissipated and reversible energies are:
\begin{align}
E_\mathrm{diss} &= \frac{1}{2} \sigma_\mathrm{c}\, \delta_\mathrm{max}\\[1ex]
E_\mathrm{rev} &= \frac{1}{2} T\, \delta
\end{align}
Moreover, a damage parameter $D$ can be defined as:
\begin{equation}
D = \min \left(
\frac{\delta_\mathrm{max}}{\delta_\mathrm{c}},1 \right)
\end{equation}
which varies from 0 (undamaged condition) and 1 (fully
damaged condition). This variable can only increase because damage is
an irreversible process. A simple penalty contact model has been incorporated
in the cohesive law so that normal tractions can be returned in
case of compression:
\begin{equation}
T_\mathrm{n} = \alpha \Delta_\mathrm{n} \quad\text{if
$\Delta_\mathrm{n} < 0$}
\end{equation}
where $\alpha$ is a stiffness parameter that defaults to zero. The
relative contact energy is equivalent to reversible energy but in
compression.
The material name of the linear decreasing cohesive law is
\code{material\_cohesive\_linear} and its parameters with their
respective default values are:
\begin{itemize}
\item \code{sigma\_c}: 0
\item \code{beta}: 0
\item \code{G\_cI}: 0
\item \code{G\_cII}: 0
\item \code{kappa}: 1
\item \code{penalty}: 0
\end{itemize}
A random number generator can be used to assign a random
$\sigma_\mathrm{c}$ to each facet following a given
distribution (see Section~\ref{sect:smm:CL}).
The bilinear constitutive law works exactly the same way as the linear
one, except for the additional parameter \code{delta\_0} that by
default is zero. Two examples for the extrinsic and intrinsic cohesive
elements and also an example to assign different properties to
intergranular and transgranular cohesive elements can be found in
the folder \code{\examplesdir/cohesive\_element/}.
-\subsubsection{Exponential Cohesive Law}
+\subsubsection{Exponential Cohesive Law}\label{ssect:smm:cl:coh-exponential}
Ortiz and Pandolfi proposed this cohesive law in 1999~\cite{ortiz1999}. The
traction-opening equation for this law is as follows:
\begin{equation}
\label{eq:exponential_law}
t = e \sigma_c \frac{\delta}{\delta_c}e^{-\delta/ \delta_c}
\end{equation}
This equation is plotted in Figure~\ref{fig:smm:CL:ECL}. The term
$\frac{\partial{(t/\delta)}}{\partial{\delta}}$ of
\eqref{eq:tangent_cohesive} for static analyses is defined as:
\begin{equation}
\frac{\partial{(t/ \delta)}}{\partial{\delta}} = \left\{\begin{array} {l l}
-e \frac{\sigma_c}{\delta_c^2 }e^{-\delta / \delta_c} & \quad if
\delta \geq \delta_{max}\\
0 & \quad if \delta < \delta_{max}, \delta_n > 0
\end{array} \right.
\end{equation}
\begin{figure}[!htb]
\begin{center}
\includegraphics[width=0.6\textwidth,keepaspectratio=true]{figures/cohesive_exponential.pdf}
\caption{Exponential cohesive law}
\label{fig:smm:CL:ECL}
\end{center}
\end{figure}
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "manual"
%%% End:
diff --git a/doc/manual/manual-constitutive-laws.tex b/doc/manual/manual-constitutive-laws.tex
index b0fb94362..087b20a5e 100644
--- a/doc/manual/manual-constitutive-laws.tex
+++ b/doc/manual/manual-constitutive-laws.tex
@@ -1,255 +1,378 @@
-\subsection{Neo-Hookean}\index{Material!Neohookean}
+\section{Constitutive Laws \label{sect:smm:CL}}\index{Material}
+In order to compute an element's response to deformation, one needs to
+use an appropriate constitutive relationship. The constitutive law is
+used to compute the element's stresses from the element's strains.
+
+In the finite-element discretization, the constitutive formulation is
+applied to every quadrature point of each element. When the implicit
+formulation is used, the tangent matrix has to be computed.
+
+The chosen materials for the simulation have to be specified in the
+mesh file or, as an alternative, they can be assigned using the
+\code{element\_material} vector. For every material assigned to the
+problem one has to specify the material characteristics (constitutive
+behavior and material properties) in a text file (\eg material.dat) as
+follows:
+\begin{cpp}
+ material %\emph{constitutive\_law}% %\emph{<optional flavor>}% [
+ name = $value$
+ rho = $value$
+ ...
+ ]
+\end{cpp}
+\index{Constitutive\_laws} where \emph{constitutive\_law} is the adopted
+constitutive law, followed by the material properties listed one by line in the
+bracket (\eg \code{name} and density \code{rho}). Some constitutive laws can
+also have an \emph{optional flavor}. For example a non-local constitutive law
+can be flavored by a weight function. The file needs to be loaded in \akantu
+using the \code{initialize} method of \akantu (as shown in
+Section~\ref{sec:writing_main})
+\begin{cpp}
+ initialize("material.dat", argc, argv);
+\end{cpp}
+% or, alternatively, the \code{initFull} method.
+% \begin{cpp}
+% model.initFull("material.dat");
+% \end{cpp}
+
+In order to conveniently store values at each quadrature in a material
+point \akantu provides a special data structure, the
+\code{InternalField}. The internal fields are inheriting from the
+\code{ElementTypeMapArray}. Furthermore, it provides several functions for
+initialization, auto-resizing and auto removal of quadrature points.
+
+Sometimes it is also desired to generate random distributions of
+internal parameters. An example might be the critical stress at which the
+material fails. To generate such a field, in the material input file,
+a random quantity needs be added to the base value:
+\begin{cpp}
+ sigma_c = $base$
+ sigma_c = $base$ uniform [$min$, $max$]
+ sigma_c = $base$ weibull [$\lambda$, $m$]
+\end{cpp}
+
+All parameters are real numbers. For the uniform distribution, minimum
+and maximum values have to be specified.
+Random parameters are defined as a $base$ value to which we add a random number
+that follows the chosen distribution.
+
+The
+\href{http://en.wikipedia.org/wiki/Uniform\_distribution\_(continuous)}{\emph{Uniform}}
+distribution is gives a random values between in $[min, max)$
+
+The \href{http://en.wikipedia.org/wiki/Weibull\_distribution}{\emph{Weibull}}
+distribution is characterized by the following cumulative distribution function:
+\begin{equation}
+ F(x) = 1- e^{-\left({x/\lambda}\right)^m}
+\end{equation}
+which depends on $m$ and $\lambda$, which are the shape parameter and the scale
+parameter.
+
+The following sections describe the constitutive models implemented in
+\akantu. In Appendix~\ref{app:material-parameters} a summary of the
+parameters for all materials of \akantu is provided.
+
+
+\subsection{Elasticity\label{ssect:smm:linear-elastic}}\index{Material!Elastic}
+The elastic law is a commonly used constitutive relationship that can
+be used for a wide range of engineering materials (\eg metals,
+concrete, rock, wood, glass, rubber, etc.) provided that the strains
+remain small (\ie small deformation and stress lower than yield
+strength). The linear elastic behavior is described by Hooke's law,
+which states that the stress is linearly proportional to the applied
+strain (material behaves like an ideal spring), as illustrated in
+Figure~\ref{fig:smm:cl:elastic}.
+\begin{figure}[!htb]
+ \begin{center}
+
+ \subfloat[]{
+ \includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/stress_strain_el.pdf}
+ \label{fig:smm:cl:elastic:stress_strain} }
+ \hspace{0.05\textwidth} \subfloat[]{
+ \raisebox{0.125\textwidth}{\includegraphics[width=0.25\textwidth,keepaspectratio=true]{figures/hooke_law.pdf}}
+ \label{fig:smm:cl:elastic:hooke} }
+ \caption{(a) Stress-strain curve for elastic material and (b)
+ schematic representation of Hooke's law, denoted as a spring.}
+ \label{fig:smm:cl:elastic}
+ \end{center}
+\end{figure}
+The equation that relates the strains to the
+displacements is: % First the strain is computed (at every gauss
+point) from the displacements as follows:
+\begin{equation}
+ \label{eqn:smm:strain_inf}
+ \mat{\epsilon} =
+ \frac{1}{2} \left[ \mat{F}-\mat{I}+(\mat{F}-\mat{I})^T \right]
+\end{equation}
+where $\mat{\epsilon}$ represents the infinitesimal strain tensor,
+$\mat{F} = \nabla_{\!\!\vec{X}}\vec{x}$ the deformation gradient
+tensor and $\mat{I}$ the identity matrix. The constitutive equation
+for isotropic homogeneous media can be expressed as:
+\begin{equation}
+ \label{eqn:smm:material:constitutive_elastic}
+ \mat{\sigma } =\lambda\mathrm{tr}(\mat{\epsilon})\mat{I}+2 \mu\mat{\epsilon}
+\end{equation}
+where $\mat{\sigma}$ is the Cauchy stress tensor
+($\lambda$ and $\mu$ are the the first and second Lame's
+coefficients). Besides the name and density, one has to specify the
+following properties in the material file: \code{E} (Young's modulus),
+\code{nu} (Poisson's ratio) and (for 2D) \code{Plane\_Stress} (if set
+to zero or not specified plane strain conditions are assumed for a
+plain analysis, otherwise plane stress conditions are applied).
+
+\subsection{Neo-Hookean\label{ssect:smm:cl:neohookean}}\index{Material!Neohookean}
The hyperelastic Neo-Hookean constitutive law results from an
extension of the linear elastic relationship (Hooke's Law) for large
deformation. Thus, the model predicts nonlinear stress-strain behavior
for bodies undergoing large deformations.
\begin{figure}[!htb]
\begin{center}
\includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/stress_strain_neo.pdf}
\caption{Neo-hookean Stress-strain curve.}
\label{fig:smm:cl:neo_hookean}
\end{center}
\end{figure}
As illustrated in Figure~\ref{fig:smm:cl:neo_hookean}, the behavior is initially
linear and the mechanical behavior is very close to the corresponding linear
elastic material. This constitutive relationship, which accounts for compressibility,
is a modified version of the one proposed by Ronald Rivlin \cite{Belytschko:2000}.
The strain energy stored in the material is given by:
\begin{equation}\label{eqn:smm:constitutive:neohookean_potential}
\Psi(\mat{C}) = \frac{1}{2}\lambda_0\left(\ln J\right)^2-\mu_0\ln J+\frac{1}{2}
\mu_0\left(\mathrm{tr}(\mat{C})-3\right)
\end{equation}
\noindent where $\lambda_0$ and $\mu_0$ are, respectively, Lamé's first parameter
and the shear modulus at the initial configuration. $J$ is the jacobian of the deformation
gradient ($\mat{F}=\nabla_{\!\!\vec{X}}\vec{x}$): $J=\text{det}(\mat{F})$. Finally $\mat{C}$ is the right Cauchy-Green
deformation tensor.
Since this kind of material is used for large deformation problems, a
finite deformation framework should be used. Therefore, the Cauchy
stress ($\mat{\sigma}$) should be computed through the second
Piola-Kirchhoff stress tensor $\mat{S}$:
\begin{equation}
\mat{\sigma } = \frac{1}{J}\mat{F}\mat{S}\mat{F}^T
\end{equation}
Finally the second Piola-Kirchhoff stress tensor is given by:
\begin{equation}
\mat{S} = 2\frac{\partial\Psi}{\partial\mat{C}} = \lambda_0\ln J
\mat{C}^{-1}+\mu_0\left(\mat{I}-\mat{C}^{-1}\right)
\end{equation}
The parameters to indicate in the material file are the same
as those for the elastic case: \code{E} (Young's modulus), \code{nu} (Poisson's
ratio).
-\subsection{Small-Deformation Plasticity}\index{Material!Small-deformation Plasticity}
+\subsection{Visco-Elasticity}
+
+\label{ssect:smm:cl:sls}
+% Standard Solid rheological model, see [] J.C. Simo, T.J.R. Hughes,
+% "Computational Inelasticity", Springer (1998), see Sections 10.2 and 10.3
+Visco-elasticity is characterized by strain rate dependent
+behavior. Moreover, when such a material undergoes a deformation it
+dissipates energy. This dissipation results in a hysteresis loop in
+the stress-strain curve at every loading cycle (see
+Figure~\ref{fig:smm:cl:visco-elastic:hyst}). In principle, it can be
+applied to many materials, since all materials exhibit a visco-elastic
+behavior if subjected to particular conditions (such as high
+temperatures).
+\begin{figure}[!htb]
+ \begin{center}
+ \subfloat[]{
+ \includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/stress_strain_visco.pdf}
+ \label{fig:smm:cl:visco-elastic:hyst}
+ }
+ \hspace{0.05\textwidth}
+ \subfloat[]{
+ \raisebox{0.025\textwidth}{\includegraphics[width=0.3\textwidth,keepaspectratio=true]{figures/visco_elastic_law.pdf}}
+ \label{fig:smm:cl:visco-elastic:model}
+ }
+ \caption{(a) Characteristic stress-strain behavior of a visco-elastic material with hysteresis loop and (b) schematic representation of the standard rheological linear solid visco-elastic model.}
+ \label{fig:smm:cl:visco-elastic}
+ \end{center}
+\end{figure}
+The standard rheological linear solid model (see Sections 10.2 and 10.3
+of~\cite{simo92}) has been implemented in \akantu. This model results from the
+combination of a spring mounted in parallel with a spring and a dashpot
+connected in series, as illustrated in
+Figure~\ref{fig:smm:cl:visco-elastic:model}. The advantage of this model is that
+it allows to account for creep or stress relaxation. The equation that relates
+the stress to the strain is (in 1D):
+\begin{equation}
+ \frac{d\epsilon(t)}{dt} = \left ( E + E_V \right ) ^ {-1} \cdot \left [ \frac{d\sigma(t)}{dt} + \frac{E_V}{\eta}\sigma(t) - \frac{EE_V}{\eta}\epsilon(t) \right ]
+\end{equation}
+where $\eta$ is the viscosity. The equilibrium condition is unique and
+is attained in the limit, as $t \to \infty $. At this stage, the
+response is elastic and depends on the Young's modulus $E$. The
+mandatory parameters for the material file are the following:
+\code{rho} (density), \code{E} (Young's modulus), \code{nu} (Poisson's
+ratio), \code{Plane\_Stress} (if set to zero plane strain, otherwise
+plane stress), \code{eta} (dashpot viscosity) and \code{Ev} (stiffness
+of the viscous element).
+
+Note that the current standard linear solid model is applied only on the deviatoric part of the strain tensor. The spheric part of the strain tensor affects the stress tensor like an linear elastic material.
+
+\subsection{Small-Deformation Plasticity\label{ssect:smm:cl:plastic}}\index{Material!Small-deformation Plasticity}
The small-deformation plasticity is a simple plasticity material
formulation which accounts for the additive decomposition of strain
into elastic and plastic strain components. This formulation is
applicable to infinitesimal deformation where the additive
decomposition of the strain is a valid approximation. In this
formulation, plastic strain is a shearing process where hydrostatic
stress has no contribution to plasticity and consequently plasticity
does not lead to volume change. Figure~\ref{fig:smm:cl:Lin-strain-hard}
shows the linear strain hardening elasto-plastic behavior according to
the additive decomposition of strain into the elastic and plastic
parts in infinitesimal deformation as
\begin{align}
\mat{\varepsilon} &= \mat{\varepsilon}^e +\mat{\varepsilon}^p\\
{\mat{\sigma}} &= 2G(\mat{\varepsilon}^e) + \lambda \mathrm{tr}(\mat{\varepsilon}^e)\mat{I}
\end{align}
\begin{figure}[htp]
\centering
{\includegraphics[scale=0.4, clip]{figures/isotropic_hardening_plasticity.pdf}}
\caption{
Stress-strain curve for the small-deformation plasticity with linear isotropic hardening.
}
\label{fig:smm:cl:Lin-strain-hard}
\end{figure}
\noindent In this class, the von Mises yield criterion is used. In the von Mises yield criterion, the yield is independent of the hydrostatic stress. Other yielding criteria such as Tresca and Gurson can be easily implemented in this class as well.
In the von Mises yield criterion, the hydrostatic stresses have no effect on the plasticity and consequently the yielding occurs when a critical elastic shear energy is achieved.
\begin{equation} \label{eqn:smm:constitutive:von Mises}
f = \sigma_{\st{eff}} - \sigma_y = \left(\frac{3}{2} {\mat{\sigma}}^{\st{tr}} : {\mat{\sigma}}^{\st{tr}}\right)^\frac{1}{2}-\sigma_y (\mat{\varepsilon}^p)
\end{equation}
\begin{equation} \label{eqn:smm:constitutive:yielding}
f < 0 \quad \textrm{Elastic deformation,} \qquad f = 0 \quad \textrm{Plastic deformation}
\end{equation}
where $\sigma_y$ is the yield strength of the material which can be function of plastic strain in case of hardening type of materials and ${\mat{\sigma}}^{\st{tr}}$ is the deviatoric part of stress given by
\begin{equation} \label{eqn:smm:constitutive:deviatoric stress}
{\mat{\sigma}}^{\st{tr}}=\mat{\sigma} - \frac{1}{3} \mathrm{tr}(\mat{\sigma}) \mat {I}
\end{equation}
After yielding $(f = 0)$, the normality hypothesis of plasticity determines the direction of plastic flow which is normal to the tangent to the yielding surface at the load point. Then, the tensorial form of the plastic constitutive equation using the von Mises yielding criterion (see equation 4.34) may be written as
\begin{equation} \label{eqn:smm:constitutive:plastic contitutive equation}
\Delta {\mat{\varepsilon}}^p = \Delta p \frac {\partial{f}}{\partial{\mat \sigma}}=\frac{3}{2} \Delta p \frac{{\mat{\sigma}}^{\st{tr}}}{\sigma_{\st{eff}}}
\end{equation}
In these expressions, the direction of the plastic strain increment (or equivalently, plastic strain rate) is given by $\frac{{\mat{\sigma}}^{\st{tr}}}{\sigma_{\st{eff}}}$ while the magnitude is defined by the plastic multiplier $\Delta p$. This can be obtained using the \emph{consistency condition} which impose the requirement for the load point to remain on the yielding surface in the plastic regime.
Here, we summarize the implementation procedures for the
small-deformation plasticity with linear isotropic hardening:
\begin{enumerate}
\item Compute the trial stress:
\begin{equation}
{\mat{\sigma}}^{\st{tr}} = {\mat{\sigma}}_t + 2G\Delta \mat{\varepsilon} + \lambda \mathrm{tr}(\Delta \mat{\varepsilon})\mat{I}
\end{equation}
\item Check the Yielding criteria:
\begin{equation}
f = (\frac{3}{2} {\mat{\sigma}}^{\st{tr}} : {\mat{\sigma}}^{\st{tr}})^{1/2}-\sigma_y (\mat{\varepsilon}^p)
\end{equation}
\item Compute the Plastic multiplier:
\begin{align}
d \Delta p &= \frac{\sigma^{tr}_{eff} - 3G \Delta P^{(k)}- \sigma_y^{(k)}}{3G + h}\\
\Delta p^{(k+1)} &= \Delta p^{(k)}+ d\Delta p\\
\sigma_y^{(k+1)} &= (\sigma_y)_t+ h\Delta p
\end{align}
\item Compute the plastic strain increment:
\begin{equation}
\Delta {\mat{\varepsilon}}^p = \frac{3}{2} \Delta p \frac{{\mat{\sigma}}^{\st{tr}}}{\sigma_{\st{eff}}}
\end{equation}
\item Compute the stress increment:
\begin{equation}
{\Delta \mat{\sigma}} = 2G(\Delta \mat{\varepsilon}-\Delta \mat{\varepsilon}^p) + \lambda \mathrm{tr}(\Delta \mat{\varepsilon}-\Delta \mat{\varepsilon}^p)\mat{I}
\end{equation}
\item Update the variables:
\begin{align}
{\mat{\varepsilon^p}} &= {\mat{\varepsilon}}^p_t+{\Delta {\mat{\varepsilon}}^p}\\
{\mat{\sigma}} &= {\mat{\sigma}}_t+{\Delta \mat{\sigma}}
\end{align}
\end{enumerate}
We use an implicit integration technique called \emph{the radial
return method} to obtain the plastic multiplier. This method has the
advantage of being unconditionally stable, however, the accuracy
remains dependent on the step size. The plastic parameters to indicate
in the material file are: \code{$\sigma_y$} (Yield stress) and
\code{h} (Hardening modulus). In addition, the elastic parameters need
to be defined as previously mentioned: \code{E} (Young's modulus),
\code{nu} (Poisson's ratio).
-
-\subsection{Visco-Elasticity}
-
-% Standard Solid rheological model, see [] J.C. Simo, T.J.R. Hughes,
-% "Computational Inelasticity", Springer (1998), see Sections 10.2 and 10.3
-Visco-elasticity is characterized by strain rate dependent
-behavior. Moreover, when such a material undergoes a deformation it
-dissipates energy. This dissipation results in a hysteresis loop in
-the stress-strain curve at every loading cycle (see
-Figure~\ref{fig:smm:cl:visco-elastic:hyst}). In principle, it can be
-applied to many materials, since all materials exhibit a visco-elastic
-behavior if subjected to particular conditions (such as high
-temperatures).
-\begin{figure}[!htb]
- \begin{center}
-
- \subfloat[]{
- \includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/stress_strain_visco.pdf}
- \label{fig:smm:cl:visco-elastic:hyst}
- }
- \hspace{0.05\textwidth}
- \subfloat[]{
- \raisebox{0.025\textwidth}{\includegraphics[width=0.3\textwidth,keepaspectratio=true]{figures/visco_elastic_law.pdf}}
- \label{fig:smm:cl:visco-elastic:model}
- }
- \caption{(a) Characteristic stress-strain behavior of a visco-elastic material with hysteresis loop and (b) schematic representation of the standard rheological linear solid visco-elastic model.}
- \label{fig:smm:cl:visco-elastic}
- \end{center}
-\end{figure}
-The standard rheological linear solid model (see Sections 10.2 and 10.3
-of~\cite{simo92}) has been implemented in \akantu. This model results from the
-combination of a spring mounted in parallel with a spring and a dashpot
-connected in series, as illustrated in
-Figure~\ref{fig:smm:cl:visco-elastic:model}. The advantage of this model is that
-it allows to account for creep or stress relaxation. The equation that relates
-the stress to the strain is (in 1D):
-\begin{equation}
- \frac{d\epsilon(t)}{dt} = \left ( E + E_V \right ) ^ {-1} \cdot \left [ \frac{d\sigma(t)}{dt} + \frac{E_V}{\eta}\sigma(t) - \frac{EE_V}{\eta}\epsilon(t) \right ]
-\end{equation}
-where $\eta$ is the viscosity. The equilibrium condition is unique and
-is attained in the limit, as $t \to \infty $. At this stage, the
-response is elastic and depends on the Young's modulus $E$. The
-mandatory parameters for the material file are the following:
-\code{rho} (density), \code{E} (Young's modulus), \code{nu} (Poisson's
-ratio), \code{Plane\_Stress} (if set to zero plane strain, otherwise
-plane stress), \code{eta} (dashpot viscosity) and \code{Ev} (stiffness
-of the viscous element).
-
-Note that the current standard linear solid model is applied only on the deviatoric part of the strain tensor. The spheric part of the strain tensor affects the stress tensor like an linear elastic material.
-
\subsection{Damage}
In the simplified case of a linear elastic and brittle material, isotropic
damage can be represented by a scalar variable $d$, which varies from $0$ to $1$
for no damage to fully broken material respectively. The stress-strain
relationship then becomes:
\begin{equation*}
\mat{\sigma} = (1-d)\, \mat{C}:\mat{\varepsilon}
\end{equation*}
where $\mat{\sigma}$, $\mat{\varepsilon}$ are the Cauchy stress and strain
tensors, and $\mat{C}$ is the elastic stiffness tensor. This formulation relies
on the definition of an evolution law for the damage variable. In \akantu, many
possibilities exist and they are listed below.
\subsubsection{Marigo}
+\label{ssect:smm:cl:damage-marigo}
This damage evolution law is energy based as defined by Marigo \cite{marigo81a,
lemaitre96a}. It is an isotropic damage law.
\begin{align}
Y &= \frac{1}{2}\mat{\varepsilon}:\mat{C}:\mat{\varepsilon}\\
F &= Y - Y_d - S d\\
d &= \left\{
\begin{array}{l l}
\mathrm{min}\left(\frac{Y-Y_d}{S},\;1\right) & \mathrm{if}\; F > 0\\
\mathrm{unchanged} & \mathrm{otherwise}
\end{array}
\right.
\end{align}
In this formulation, $Y$ is the strain energy release rate, $Y_d$ the
rupture criterion and $S$ the damage energy. The non-local version of
this damage evolution law is constructed by averaging the energy $Y$.
\subsubsection{Mazars}
+\label{ssect:smm:cl:damage-mazars}
This law introduced by Mazars \cite{mazars84a} is a behavioral model to
represent damage evolution in concrete. The governing variable in this damage
law is the equivalent strain $\varepsilon_{\st{eq}} =
\sqrt{<\mat{\varepsilon}>_+:<\mat{\varepsilon}>_+}$, with $<.>_+$ the positive
part of the tensor.
The damage the is defined as:
\begin{align}
D &= \alpha_t^\beta D_t + (1-\alpha_t)^\beta D_c\\
D_t &= 1 - \frac{\kappa_0 (1- A_t)}{\varepsilon_{\st{eq}}} - A_t \exp^{-B_t(\varepsilon_{\st{eq}}-\kappa_0)}\\
D_c &= 1 - \frac{\kappa_0 (1- A_c)}{\varepsilon_{\st{eq}}} - A_c
\exp^{-B_c(\varepsilon_{\st{eq}}-\kappa_0)}\\
\alpha_t &= \frac{\sum_{i=1}^3<\varepsilon_i>_+\varepsilon_{\st{nd}\;i}}{\varepsilon_{\st{eq}}^2}
\end{align}
With $\kappa_0$ the damage threshold, $A_t$ and $B_t$ the damage parameter in
traction, $A_c$ and $B_c$ the damage parameter in compression, $\beta$ is the
shear parameter. $\alpha_t$ is the coupling parameter between traction and
compression, the $\varepsilon_i$ are the eigenstrain and the
$\varepsilon_{\st{nd}\;i}$ are the eigenvalues of the strain if the material
were undamaged.
The coefficients $A$ and $B$ are the post-peak asymptotic
value and the decay shape parameters.
\IfFileExists{manual-extra_materials.tex}{\input{manual-extra_materials}}{}
\IfFileExists{manual-cohesive_laws.tex}{\input{manual-cohesive_laws}}{}
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "manual"
%%% End:
diff --git a/doc/manual/manual-solidmechanicsmodel.tex b/doc/manual/manual-solidmechanicsmodel.tex
index 553b6f6e8..67c6c7ef9 100644
--- a/doc/manual/manual-solidmechanicsmodel.tex
+++ b/doc/manual/manual-solidmechanicsmodel.tex
@@ -1,1127 +1,1020 @@
\chapter{Solid Mechanics Model\index{SolidMechanicsModel}\label{sect:smm}}
The solid mechanics model is a specific implementation of the
\code{Model} interface dedicated to handle the equations of motion or
equations of equilibrium. The model is created for a given mesh. It
will create its own \code{FEEngine} object to compute the interpolation,
gradient, integration and assembly operations. A
\code{SolidMechanicsModel} object can simply be created like this:
\begin{cpp}
SolidMechanicsModel model(mesh, spatial_dimension);
\end{cpp}
where \code{mesh} is the mesh for which the equations are to be
solved, and \code{spatial\_dimension} is the dimensionality of the
problem. If the \code{spatial\_dimension} is omitted, the problem is
assumed to have the same dimensionality as the one specified in the
mesh.
This model contains at least the the following six \code{Arrays}:
\begin{description}
\item[blocked\_dofs] contains a Boolean value for each degree of
freedom specifying whether that degree is blocked or not. A
Dirichlet boundary condition can be prescribed by setting the
\textbf{blocked\_dofs} value of a degree of freedom to \code{true}.
A Neumann boundary condition can be applied by setting the
\textbf{blocked\_dofs} value of a degree of freedom to \code{false}.
The \textbf{displacement}, \textbf{velocity} and
\textbf{acceleration} are computed for all degrees of freedom for
which the \textbf{blocked\_dofs} value is set to \code{false}. For
the remaining degrees of freedom, the imposed values (zero by
default after initialization) are kept.
\item[displacement] contains the displacements of all degrees of
freedom. It can be either a computed displacement for free degrees
of freedom or an imposed displacement in case of blocked ones
($\vec{u}$ in the following).
\item[velocity] contains the velocities of all degrees of freedom. As
\textbf{displacement}, it contains computed or imposed velocities
depending on the nature of the degrees of freedom ($\vec{\dot{u}}$
in the following).
\item[acceleration] contains the accelerations of all degrees of
freedom. As \textbf{displacement}, it contains computed or imposed
accelerations depending on the nature of the degrees of freedom
($\vec{\ddot{u}}$ in the following).
\item[force] contains the external forces applied on the nodes
($\vec{f_{\st{ext}}}$ in the following).
\item[residual] contains the difference between external and internal
forces. On blocked degrees of freedom, \textbf{residual} contains
the support reactions. ($\vec{r}$ in the following). It should be
mentioned that at equilibrium \textbf{residual} should be zero on
free degrees of freedom.
\end{description}
Some examples to help to understand how to use this model will be
presented in the next sections.
\section{Model Setup}
\subsection{Setting Initial Conditions \label{sect:smm:initial_condition}}
For a unique solution of the equations of motion, initial
displacements and velocities for all degrees of freedom must be
specified:
\begin{eqnarray}
\vec{u(t=0)} = \vec{u}_{0}\\
\vec{\dot u(t=0)} =\vec{v}_{0}
\end{eqnarray} The solid mechanics model can be initialized as
follows:
\begin{cpp}
model.initFull()
\end{cpp}
This function initializes the internal arrays and sets them to
zero. Initial displacements and velocities that are not equal to zero
can be prescribed by running a loop over the total number of
nodes. Here, the initial displacement in $x$-direction and the
initial velocity in $y$-direction for all nodes is set to $0.1$ and $1$,
respectively.
\begin{cpp}
Array<Real> & disp = model.getDisplacement();
Array<Real> & velo = model.getVelocity();
for (UInt i = 0; i < nb_nodes; ++i) {
disp(i, 0) = 0.1;
velo(i, 1) = 1.;
}
\end{cpp}
\subsection{Setting Boundary Conditions\label{sect:smm:boundary}}
This section explains how to impose Dirichlet or Neumann boundary
conditions. A Dirichlet boundary condition specifies the values that
the displacement needs to take for every point $x$ at the boundary
($\Gamma_u$) of the problem domain (Fig.~\ref{fig:smm:boundaries}):
\begin{equation}
\vec{u} = \vec{\bar u} \quad \forall \vec{x}\in
\Gamma_{u}
\end{equation}
A Neumann boundary condition imposes the value of the gradient of the
solution at the boundary ($\Gamma_t$) of the problem domain
(Fig.~\ref{fig:smm:boundaries}):
\begin{equation}
\vec{t} = \mat{\sigma} \vec{n} = \vec{\bar t} \quad
\forall \vec{x}\in \Gamma_{t}
\end{equation}
\begin{figure} \centering
\def\svgwidth{0.5\columnwidth}
\input{figures/problemDomain.pdf_tex}
\caption{Problem domain $\Omega$ with boundary in three dimensions. The Dirchelet and the Neumann regions of the boundary are denoted with $\Gamma_u$ and $\Gamma_t$, respecitvely.\label{fig:smm:boundaries}}
\label{fig:problemDomain}
\end{figure}
Different ways of imposing these boundary conditions exist. A basic
way is to loop over nodes or elements at the boundary and apply local
values. A more advanced method consists of using the notion of the
boundary of the mesh. In the following both ways are presented.
Starting with the basic approach, as mentioned, the Dirichlet boundary
conditions can be applied by looping over the nodes and assigning the
required values. Figure~\ref{fig:smm:dirichlet_bc} shows a beam with a
fixed support on the left side. On the right end of the beam, a load
is applied. At the fixed support, the displacement has a given
value. For this example, the displacements in both the $x$ and the
$y$-direction are set to zero. Implementing this displacement boundary
condition is similar to the implementation of initial displacement
conditions described above. However, in order to impose a displacement
boundary condition for all time steps, the corresponding nodes need to
be marked as boundary nodes as shown in the following code:
\begin{cpp}
Array<bool> & blocked = model.getBlockedDOFs();
const Array<Real> & pos = mesh.getNodes();
UInt nb_nodes = mesh.getNbNodes();
Real epsilon = Math::getTolerance();
/* this tolerance is by default equal to the machine precision but can be changed by using %\code{Math::setTolerance(value)}% */
for (UInt i = 0; i < nb_nodes; ++i) {
if(std::abs(pos(i, 0)) < epsilon) {
blocked(i, 0) = true; //block displacement in x-direction
blocked(i, 1) = true; //block displacement in y-direction
disp(i, 0) = 0.; //fixed displacement in x-direction
disp(i, 1)= 0.; //fixed displacement in y-direction
}
}
\end{cpp}
\begin{figure}[!htb]
\centering
\includegraphics[scale=0.4]{figures/dirichlet}
\caption{Beam with fixed support.\label{fig:smm:dirichlet_bc}}
\end{figure}
For the more advanced approach, one needs the notion of a boundary in
the mesh. Therefore, the boundary should be created before boundary
condition functors can be applied. Generally the boundary can be
specified from the mesh file or the geometry. For the first case, the
function \code{createGroupsFromMeshData} is called. This function
can read any types of mesh data which are provided in the mesh
file. If the mesh file is created with Gmsh, the function takes one
input strings which is either \code{tag\_0}, \code{tag\_1} or
\code{physical\_names}. The first two tags are assigned by Gmsh to
each element which shows the physical group that they belong to. In
Gmsh, it is also possible to consider strings for different groups of
elements. These elements can be separated by giving a string
\code{physical\_names} to the function
\code{createGroupsFromMeshData}. Boundary conditions can also be
created from the geometry by calling
\code{createBoundaryGroupFromGeometry}. This function gathers all the
elements on the boundary of the geometry.
To apply the required boundary conditions, the function \code{applyBC}
from the \code{SolidMechanicsModel} needs to be called. This function
gets a Dirichlet or Neumann functor and a string which specifies the
desired boundary on which the boundary conditions is to be
applied. The functors specify the type of conditions to apply. Three
built-in functors for Dirichlet exist: \code{FlagOnly, FixedValue,}
and \code{IncrementValue}. The functor \code{FlagOnly} is used if a
point is fixed in a given direction. Therefore, the input parameter to
this functor is only the fixed direction. The \code{FixedValue}
functor is used when a displacement value is applied in a fixed
direction. The \code{IncrementValue} applies an increment to the
displacement in a given direction. The following code shows the
utilization of three functors for the top, bottom and side surface of
the mesh which were already defined in the Gmsh file:
\begin{cpp}
mesh.createGroupsFromMeshData<std::string>("physical_names");
model.applyBC(BC::Dirichlet::FixedValue(13.0, BC::_y), "Top");
model.applyBC(BC::Dirichlet::FlagOnly(BC::_x), "Bottom");
model.applyBC(BC::Dirichlet::IncrementValue(13.0, BC::_x), "Side");
\end{cpp}
To apply a Neumann boundary condition, the applied traction or stress
should be specified before. In case of specifying the traction on the
surface, the functor \code{FromTraction} of Neumann boundary
conditions is called. Otherwise, the functor \code{FromStress} should
be called which gets the stress tensor as an input parameter.
\begin{cpp}
Array<Real> surface_traction(3);
surface_traction(0)=0.0;
surface_traction(1)=0.0;
surface_traction(2)=-1.0;
Array<Real> surface_stress(3, 3, 0.0);
surface_stress(0,0)=0.0;
surface_stress(1,1)=0.0;
surface_stress(2,2)=-1.0;
model.applyBC(BC::Neumann::FromTraction(surface_traction), "Bottom");
model.applyBC(BC::Neumann::FromStress(surface_stress), "Top");
\end{cpp}
If the boundary conditions need to be removed during the simulation, a
functor is called from the Neumann boundary condition to free those
boundary conditions from the desired boundary.
\begin{cpp}
model.applyBC(BC::Neumann::FreeBoundary(), "Side");
\end{cpp}
User specified functors can also be implemented. A full example for
setting both initial and boundary conditions can be found in
\shellcode{\examplesdir/boundary\_conditions.cc}. The problem solved
in this example is shown in Fig.~\ref{fig:smm:bc_and_ic}. It consists
of a plate that is fixed with movable supports on the left and bottom
side. On the right side, a traction, which increases linearly with the
number of time steps, is applied. The initial displacement and
velocity in $x$-direction at all free nodes is zero and two
respectively.
\begin{figure}[!htb]
\centering
\includegraphics[scale=0.8]{figures/bc_and_ic_example}
\caption{Plate on movable supports.\label{fig:smm:bc_and_ic}}
\end{figure}
\subsection{Material Selector\label{sect:smm:materialselector}}
If the user wants to assign different materials to the different
finite elements of choice in \akantu, a material selector has to be
used. There are four different material selectors pre-defined in
\akantu.
\begin{enumerate}
\item \code{MaterialSelector}: This material selector assigns a
material to a specific element using an index value. If an index is
not specified, it is using \code{0}, hence the first material as
default.
\item \code{DefaultMaterialSelector}: This selector assigns materials
based on the values in the \code{element\_index\_by\_material}
array.
\item \code{MeshDataMaterialSelector}: This material selector class
uses mesh data information to assign different materials. With the
proper tag name and index, different materials can be assigned as
demonstrated in the example below.
\item \code{DefaultMaterialCohesiveSelector}: This is the default
material selector for the cohesive elements and it inherits its
properties from \code{DefaultMaterialSelector}. This material
selector first checks if a cohesive element is present on an element
facet. If there is, then the first cohesive element in the material
file is assigned to that element. Otherwise, this material is
assigned to the element facet, to which, a cohesive element might be
inserted later. Hence, this material selector works for the
extrinsic cohesive elements which are dynamically inserted
throughout the analysis. If the element of interest is not a
cohesive element, then \code{DefaultMaterialSelector} is used by
default.
\end{enumerate}
Apart from the \akantu's default material selectors, users can always
develop their own classes in the main code to tackle various
multi-material assignment situations.
For example, an application of \code{MeshDataMaterialSelector} could
look like this:
\begin{cpp}
MeshDataMaterialSelector<std::string> * mat_selector;
mat_selector = new MeshDataMaterialSelector<std::string>("physical_names", model);
model.setMaterialSelector(*mat_selector);
\end{cpp}
In this example the physical names specified in a GMSH geometry file will by
used to match the material names in the input file.
Another example would be:
\begin{cpp}
MeshDataMaterialSelector<UInt> * mat_selector;
mat_selector = new MeshDataMaterialSelector<UInt>("tag_1", model);
model.setMaterialSelector(*mat_selector);
\end{cpp}
where \code{tag\_1} of the mesh file is used as the classifier index
and the elements that have index value equal to one will be assigned
to the second material in the material file.
\IfFileExists{manual-cohesive_elements_insertion.tex}{\input{manual-cohesive_elements_insertion}}{}
\section{Static Analysis\label{sect:smm:static}}
The \code{SolidMechanicsModel} class can handle different analysis
methods, the first one being presented is the static case. In this
case, the equation to solve is,
\begin{equation}
\label{eqn:smm:static} \mat{K} \vec{u} =
\vec{f_{\st{ext}}}~,
\end{equation}
where $\mat{K}$ is the global stiffness matrix, $\vec{u}$ the
displacement vector and $\vec{f_{\st{ext}}}$ the vector of external
forces applied to the system.
To solve such a problem, the static solver of the
\code{SolidMechanicsModel}\index{SolidMechanicsModel} object is used.
First, a model has to be created and initialized. To create the
model, a mesh (which can be read from a file) is needed, as explained
in Section~\ref{sect:common:mesh}. Once an instance of a
\code{SolidMechanicsModel} is obtained, the easiest way to initialize
it is to use the \code{initFull}\index{SolidMechanicsModel!initFull}
method by giving the \code{SolidMechanicsModelOptions}. These options
specify the type of analysis to be performed and whether the materials
should be initialized or not.
\begin{cpp}
SolidMechanicsModel model(mesh);
model.initFull(SolidMechanicsModelOptions(_static));
\end{cpp}
Here, a static analysis is chosen by passing the argument
\code{\_static} to the method. By default, the Boolean for no
initialization of the materials is set to false, so that they are
initialized during the \code{initFull}. The method \code{initFull}
also initializes all appropriate vectors to zero. Once the model is
created and initialized, the boundary conditions can be set as
explained in Section~\ref{sect:smm:boundary}. Boundary conditions
will prescribe the external forces for some free degrees of freedom
$\vec{f_{\st{ext}}}$ and displacements for some others. At this point
of the analysis, the function
\code{solveStep}\index{SolidMechanicsModel!solveStep} can be called:
\begin{cpp}
model.solveStep<_scm_newton_raphson_tangent_modified, _scc_residual>(1e-4, 1);
\end{cpp}
This function is templated by the solving method and the convergence
criterion and takes two arguments: the tolerance and the maximum
number of iterations, which are $1e-4$ and $1$ for this example. The
modified Newton-Raphson method is chosen to solve the system. In this
method, the equilibrium equation (\ref{eqn:smm:static}) is modified in
order to apply a Newton-Raphson convergence algorithm:
\begin{align}\label{eqn:smm:static-newton-raphson}
\mat{K}^{i+1}\delta\vec{u}^{i+1} &= \vec{r} \\
&= \vec{f_{\st{ext}}} -\vec{f_{\st{int}}}\\
&= \vec{f_{\st{ext}}} - \mat{K}^{i} \vec{u}^{i}\\
\vec{u}^{i+1} &= \vec{u}^{i} + \delta\vec{u}^{i+1}~,\nonumber
\end{align}
where $\delta\vec{ u}$ is the increment of displacement to be added
from one iteration to the other, and $i$ is the Newton-Raphson
iteration counter. By invoking the \code{solveStep} method in the
first step, the global stiffness matrix $\mat{K}$ from
Equation~(\ref{eqn:smm:static}) is automatically assembled. A
Newton-Raphson iteration is subsequently started, $\mat{K}$ is updated
according to the displacement computed at the previous iteration and
one loops until the forces are balanced (\code{\_scc\_residual}), \ie
$||\vec{r}|| < \mbox{\code{\_scc\_residual}}$. One can also iterate
until the increment of displacement is zero (\code{\_scc\_increment})
which also means that the equilibrium is found. For a linear elastic
problem, the solution is obtained in one iteration and therefore the
maximum number of iterations can be set to one. But for a non-linear
case, one needs to iterate as long as the norm of the residual exceeds
the tolerance threshold and therefore the maximum number of iterations
has to be higher, e.g. $100$:
\begin{cpp}
model.solveStep<_scm_newton_raphson_tangent_modified,_scc_residual>(1e-4, 100);
\end{cpp}
At the end of the analysis, the final solution is stored in the
\textbf{displacement} vector. A full example of how to solve a static
problem is presented in the code \code{\examplesdir/static/static.cc}.
This example is composed of a 2D plate of steel, blocked with rollers
on the left and bottom sides as shown in Figure \ref{fig:smm:static}.
The nodes from the right side of the sample are displaced by $0.01\%$
of the length of the plate.
\begin{figure}[!htb]
\centering
\includegraphics{figures/static}
\caption{Numerical setup\label{fig:smm:static}}
\end{figure}
The results of this analysis is depicted in
Figure~\ref{fig:smm:implicit:static_solution}.
\begin{figure}[!htb]
\centering
\includegraphics[width=.6\linewidth]{figures/static_analysis}
\caption{Solution of the static analysis. Left: the initial
condition, right: the solution (deformation magnified 50 times)}
\label{fig:smm:implicit:static_solution}
\end{figure}
\section{Dynamic Methods} \label{sect:smm:Dynamic_methods}
Different ways to solve the equations of motion are implemented in the
solid mechanics model. The complete equations that should be solved
are:
\begin{equation}
\label{eqn:equation-motion}
\mat{M}\vec{\ddot{u}} +
\mat{C}\vec{\dot{u}} + \mat{K}\vec{u} = \vec{f_{\st{ext}}}~,
\end{equation}
where $\mat{M}$, $\mat{C}$ and $\mat{K}$ are the mass,
damping and stiffness matrices, respectively.
In the previous section, it has already been discussed how to solve
this equation in the static case, where $\vec{\ddot{u}} =
\vec{\dot{u}} = 0$. Here the method to solve this equation in the
general case will be presented. For this purpose, a time
discretization has to be specified. The most common discretization
method in solid mechanics is the Newmark-$\beta$ method, which is
also the default in \akantu.
For the Newmark-$\beta$ method, (\ref{eqn:equation-motion}) becomes a
system of three equations (see \cite{curnier92a} \cite{hughes-83a} for
more details):
\begin{align}
\mat{M} \vec{\ddot{u}}_{n+1} + \mat{C}\vec{\dot{u}}_{n+1} + \mat{K} \vec{u}_{n+1} &=\vec{f_{\st{ext}}}_{n+1}
\label{eqn:equation-motion-discret} \\
\vec{u}_{n+1} &=\vec{u}_{n} + \left(1 - \alpha\right) \Delta t \vec{\dot{u}}_{n} +
\alpha \Delta t \vec{\dot{u}}_{n+1} + \left(\frac{1}{2} -
\alpha\right) \Delta t^2
\vec{\ddot{u}}_{n} \label{eqn:finite-difference-1}\\
\vec{\dot{u}}_{n+1} &= \vec{\dot{u}}_{n} + \left(1 - \beta\right)
\Delta t \vec{\ddot{u}}_{n} + \beta \Delta t
\vec{\ddot{u}}_{n+1} \label{eqn:finite-difference-2}
\end{align}
In these new equations, $\vec{\ddot{u}}_{n}$, $\vec{\dot{u}}_{n}$ and
$\vec{u}_{n}$ are the approximations of $\vec{\ddot{u}(t_n)}$,
$\vec{\dot{u}(t_n)}$ and $\vec{u(t_n)}$.
Equation~(\ref{eqn:equation-motion-discret}) is the equation of motion
discretized in space (finite-element discretization), and equations
(\ref{eqn:finite-difference-1}) and (\ref{eqn:finite-difference-2})
are discretized in both space and time (Newmark discretization). The
$\alpha$ and $\beta$ parameters determine the stability and the
accuracy of the algorithm. Classical values for $\alpha$ and $\beta$
are usually $\beta = 1/2$ for no numerical damping and $0 < \alpha <
1/2$.
\begin{center}
\begin{tabular}{cll}
\toprule
$\alpha$ & Method ($\beta = 1/2$) & Type\\
\midrule
$0$ & central difference & explicit\\
$1/6$ & Fox-Goodwin (royal road) &implicit\\
$1/3$ & Linear acceleration &implicit\\
$1/2$ & Average acceleration (trapezoidal rule)& implicit\\
\bottomrule
\end{tabular}
\end{center}
The solution of this system of equations,
(\ref{eqn:equation-motion-discret})-(\ref{eqn:finite-difference-2}) is
split into a predictor and a corrector system of equations. Moreover,
in the case of a non-linear equations, an iterative algorithm such as
the Newton-Raphson method is applied. The system of equations can be
written as:
\noindent\textit{Predictor:}
\begin{align}
\vec{u}_{n+1}^{0} &= \vec{u}_{n} + \Delta t
\vec{\dot{u}}_{n} + \frac{\Delta t^2}{2} \vec{\ddot{u}}_{n} \\
\vec{\dot{u}}_{n+1}^{0} &= \vec{\dot{u}}_{n} + \Delta t
\vec{\ddot{u}}_{n} \\
\vec{\ddot{u}}_{n+1}^{0} &= \vec{\ddot{u}}_{n}
\end{align}
\noindent\textit{Solve:}
\begin{align}
\left(c \mat{M} + d \mat{C} + e \mat{K}_{n+1}^i\right)
\vec{w} = \vec{f_{\st{ext}}}_{n+1} - \vec{f_{\st{int}}}_{n+1}^i -
\mat{C} \vec{\dot{u}}_{n+1}^i - \mat{M} \vec{\ddot{u}}_{n+1}^i = \vec{r}_{n+1}^i
\end{align}
\noindent\textit{Corrector:}
\begin{align}
\vec{\ddot{u}}_{n+1}^{i+1} &= \vec{\ddot{u}}_{n+1}^{i} +c \vec{w} \\
\vec{\dot{u}}_{n+1}^{i+1} &= \vec{\dot{u}}_{n+1}^{i} + d\vec{w} \\
\vec{u}_{n+1}^{i+1} &= \vec{u}_{n+1}^{i} + e \vec{w}
\end{align}
where $i$ is the Newton-Raphson iteration counter and $c$, $d$ and $e$
are parameters depending on the method used to solve the equations
\begin{center}
\begin{tabular}{lcccc}
\toprule
& $\vec{w}$ & $e$ & $d$ & $c$\\
\midrule
in acceleration &$ \delta\vec{\ddot{u}}$ & $\alpha \beta\Delta t^2$ &$\beta \Delta t$ &$1$\\
in velocity & $ \delta\vec{\dot{u}}$& $\frac{1}{\beta} \Delta t$ & $1$ & $\alpha\Delta t$\\
in displacement &$\delta\vec{u}$ & $ 1$ & $\frac{1}{\alpha} \Delta t$ & $\frac{1}{\alpha \beta} \Delta t^2$\\
\bottomrule
\end{tabular}
\end{center}
%\note{If you want to use the implicit solver \akantu should be compiled at least with one sparse matrix solver such as Mumps\cite{mumps}.}
\subsection{Implicit Time Integration}
To solve a problem with an implicit time integration scheme, first a
\code{SolidMechanicsModel} object has to be created and initialized.
Then the initial and boundary conditions have to be set. Everything
is similar to the example in the static case
(Section~\ref{sect:smm:static}), however, in this case the implicit
dynamic scheme is selected at the initialization of the model.
\begin{cpp}
SolidMechanicsModel model(mesh);
model.initFull(SolidMechanicsModelOptions(_implicit_dynamic));
/*Boundary conditions see Section~%\ref{sect:smm:boundary}% */
\end{cpp}
Because a dynamic simulation is conducted, an integration time step
$\Delta t$ has to be specified. In the case of implicit simulations,
\akantu implements a trapezoidal rule by default. That is to say
$\alpha = 1/2$ and $\beta = 1/2$ which is unconditionally
stable. Therefore the value of the time step can be chosen arbitrarily
within reason. \index{SolidMechanicsModel!setTimeStep}
\begin{cpp}
model.setTimeStep(time_step);
\end{cpp}
Since the system has to be solved for a given amount of time, the
method \code{solveStep()}, (which has already been used in the static
example in Section~\ref{sect:smm:static}), is called inside a time
loop:
\begin{cpp}
/// time loop
Real time = 0.;
for (UInt s = 1; time <max_time; ++s, time += time_step) {
model.solveStep<_scm_newton_raphson_tangent_modified,_scc_increment>(1e-12, 100);
}
\end{cpp}
An example of solid mechanics with an implicit time integration scheme
is presented in
\shellcode{\examplesdir/implicit/implicit\_dynamic.cc}. This example
consists of a 3D beam of
$\unit{10}{\meter}\,\times\,\unit{1}{\meter}\,\times\,\unit{1}{\meter}$
blocked on one side and is on a roller on the other side. A constant
force of \unit{5}{\kilo\newton} is applied in its middle.
Figure~\ref{fig:smm:implicit:dynamic} presents the geometry of this
case. The material used is a fictitious linear elastic material with a
density of \unit{1000}{\kilogrampercubicmetre}, a Young's Modulus of
\unit{120}{\mega\pascal} and Poisson's ratio of $0.3$. These values
were chosen to simplify the analytical solution.
An approximation of the dynamic response of the middle point of the
beam is given by:
\begin{equation}
\label{eqn:smm:implicit}
u\left(\frac{L}{2}, t\right)
= \frac{1}{\pi^4} \left(1 - cos\left(\pi^2 t\right) +
\frac{1}{81}\left(1 - cos\left(3^2 \pi^2 t\right)\right) +
\frac{1}{625}\left(1 - cos\left(5^2 \pi^2 t\right)\right)\right)
\end{equation}
\begin{figure}[!htb]
\centering
\includegraphics[scale=.6]{figures/implicit_dynamic}
\caption{Numerical setup}
\label{fig:smm:implicit:dynamic}
\end{figure}
Figure \ref{fig:smm:implicit:dynamic_solution} presents the deformed
beam at 3 different times during the simulation: time steps 0, 1000 and
2000.
\begin{figure}[!htb]
\centering
\setlength{\unitlength}{0.1\textwidth}
\begin{tikzpicture}
\node[above right] (img) at (0,0)
{\includegraphics[width=.6\linewidth]{figures/dynamic_analysis}};
\node[left] at (0pt,20pt) {$0$}; \node[left] at (0pt,60pt) {$1000$};
\node[left] at (0pt,100pt) {$2000$};
\end{tikzpicture}
\caption{Deformed beam at 3 different times (displacement are
magnified by a factor 10).}
\label{fig:smm:implicit:dynamic_solution}
\end{figure}
\subsection{Explicit Time Integration}
The explicit dynamic time integration scheme is based on the
Newmark-$\beta$ scheme with $\alpha=0$ (see equations
\ref{eqn:equation-motion-discret}-\ref{eqn:finite-difference-2}). In
\akantu, $\beta$ is defaults to $\beta=1/2$, see section
\ref{sect:smm:Dynamic_methods}.
The initialization of the simulation is similar to the static and
implicit dynamic version. The model is created from the
\code{SolidMechanicsModel} class. In the initialization, the explicit
scheme is selected using the \code{\_explicit\_lumped\_mass} constant.
\begin{cpp}
SolidMechanicsModel model(mesh);
model.initFull(SolidMechanicsModelOptions(_explicit_lumped_mass));
\end{cpp}
\index{SolidMechanicsModel!initFull}
\note{Writing \code{model.initFull()} or \code{model.initFull(SolidMechanicsModelOptions());} is
equivalent to use the \code{\_explicit\_lumped\_mass} keyword, as this
is the default case.}
The explicit time integration scheme implemented in \akantu uses a
lumped mass matrix $\mat{M}$ (reducing the computational cost). This
matrix is assembled by distributing the mass of each element onto its
nodes. The resulting $\mat{M}$ is therefore a diagonal matrix stored
in the \textbf{mass} vector of the model.
The explicit integration scheme is conditionally stable. The time step
has to be smaller than the stable time step which is obtained in
\akantu as follows:
\begin{cpp}
time_step = model.getStableTimeStep();
\end{cpp} \index{SolidMechanicsModel!StableTimeStep}
The stable time step is defined as:
\begin{equation}
\label{eqn:smm:explicit:stabletime}
\Delta t_{\st{crit}} = \Delta x \sqrt{\frac{\rho}{2 \mu + \lambda}}
\end{equation}
where $\Delta x$ is a characteristic length (\eg the inradius in the
case of linear triangle element), $\mu$ and $\lambda$ are the first
and second Lame's coefficients and $\rho$ is the density. The stable
time step corresponds to the time the fastest wave (the compressive
wave) needs to travel the characteristic length of the mesh.
However, it is recommended to impose a time step that is smaller than the stable time step, for
instance, by multiplying the stable time step by a safety factor
smaller than one.
\begin{cpp}
const Real safety_time_factor = 0.8;
Real applied_time_step = time_step * safety_time_factor;
model.setTimeStep(applied_time_step);
\end{cpp}
\index{SolidMechanicsModel!setTimeStep} The initial displacement and
velocity fields are, by default, equal to zero if not given
specifically by the user (see \ref{sect:smm:initial_condition}).
Like in implicit dynamics, a time loop is used in which the
displacement, velocity and acceleration fields are updated at each
time step. The values of these fields are obtained from the
Newmark$-\beta$ equations with $\beta=1/2$ and $\alpha=0$. In \akantu
these computations at each time step are invoked by calling the
function \code{solveStep}:
\begin{cpp}
for (UInt s = 1; (s-1)*applied_time_step < total_time; ++s) {
model.solveStep();
}
\end{cpp} \index{SolidMechanicsModel!solveStep}
The method
\code{solveStep} wraps the four following functions:
\begin{itemize}
\item \code{model.explicitPred()} allows to compute the displacement
field at $t+1$ and a part of the velocity field at $t+1$, denoted by
$\vec{\dot{u}^{\st{p}}}_{n+1}$, which will be used later in the method
\code{model.explicitCorr()}. The equations are:
\begin{align}
\vec{u}_{n+1} &= \vec{u}_{n} + \Delta t
\vec{\dot{u}}_{n} + \frac{\Delta t^2}{2} \vec{\ddot{u}}_{n}\\
\vec{\dot{u}^{\st{p}}}_{n+1} &= \vec{\dot{u}}_{n} + \Delta t
\vec{\ddot{u}}_{n}
\label{eqn:smm:explicit:onehalfvelocity}
\end{align}
\item \code{model.updateResidual()} and
\code{model.updateAcceleration()} compute the acceleration increment
$\delta \vec{\ddot{u}}$:
\begin{equation}
\left(\mat{M} + \frac{1}{2} \Delta t \mat{C}\right)
\delta \vec{\ddot{u}} = \vec{f_{\st{ext}}} - \vec{f}_{\st{int}\, n+1}
- \mat{C} \vec{\dot{u}}_{n} - \mat{M} \vec{\ddot{u}}_{n}
\end{equation}
\note{The internal force $\vec{f}_{\st{int}\, n+1}$ is computed from
the displacement $\vec{u}_{n+1}$ based on the constitutive law.}
\item \code{model.explicitCorr()} computes the velocity and
acceleration fields at $t+1$:
\begin{align}
\vec{\dot{u}}_{n+1} &= \vec{\dot{u}^{\st{p}}}_{n+1} + \frac{\Delta t}{2}
\delta \vec{\ddot{u}} \\ \vec{\ddot{u}}_{n+1} &=
\vec{\ddot{u}}_{n} + \delta \vec{\ddot{u}}
\end{align}
\end{itemize}
The use of the explicit time integration scheme is illustrated by an
example (see \code{\examplesdir/explicit/explicit\_dynamic.cc}). This
example models the propagation of a wave in a steel beam. The beam is
blocked on one side and a displacement is applied on the other side,
as shown in Figure~\ref{fig:smm:explicit}.
\begin{figure}[!htb] \centering
\includegraphics[scale=.6]{figures/explicit_dynamic}
\caption{Numerical setup \label{fig:smm:explicit}}
\end{figure}
The length and height of the beam are \unit{10}{\metre} and
\unit{1}{\metre}, respectively. The material is linear elastic,
homogeneous and isotropic (density:
\unit{7800}{\kilogrampercubicmetre}, Young's modulus:
\unit{210}{\giga\pascal} and Poisson's ratio: $0.3$). The imposed
displacement is equal to $\Delta u = \unit{0.05}{\metre}$. The
potential, kinetic and total energies are computed. The safety factor
is equal to $0.1$. The total simulated time is \unit{0.01}{\second}.
-\section{Constitutive Laws \label{sect:smm:CL}}\index{Material}
-In order to compute an element's response to deformation, one needs to
-use an appropriate constitutive relationship. The constitutive law is
-used to compute the element's stresses from the element's strains.
-
-In the finite-element discretization, the constitutive formulation is
-applied to every quadrature point of each element. When the implicit
-formulation is used, the tangent matrix has to be computed.
-
-The chosen materials for the simulation have to be specified in the
-mesh file or, as an alternative, they can be assigned using the
-\code{element\_material} vector. For every material assigned to the
-problem one has to specify the material characteristics (constitutive
-behavior and material properties) in a text file (\eg material.dat) as
-follows:
-\begin{cpp}
- material %\emph{constitutive\_law}% [
- name = %\emph{XXX}%
- rho = %\emph{XXX}%
- ...
- ]
-\end{cpp}
-\index{Constitutive\_laws} where \emph{constitutive\_law} is the
-adopted constitutive law, followed by the material properties listed
-one by line in the bracket (\eg \code{name} and density
-\code{rho}). The file needs to be loaded in \akantu using the
-\code{initialize} method of \akantu (as shown in Section~\ref{sec:writing_main})
-\begin{cpp}
-initialize("material.dat", argc, argv);
-\end{cpp}
-% or, alternatively, the \code{initFull} method.
-% \begin{cpp}
-% model.initFull("material.dat");
-% \end{cpp}
-
-In order to conveniently store stresses and strains at each quadrature
-point \akantu provides a special data structure, the
-\code{InternalField}. The internal field is inheriting from the
-\code{ElementTypeMapArray}. It has a reference to the material to
-which it belongs. Furthermore, it provides several functions for
-initialization, resizing and removal of quadrature points. These
-functions are for instance used inside the materials, if elements are
-removed from one material and inserted into another one.
-
-Sometimes it is also desired to generate random distributions of
-internal values. An example might be the critical stress at which the
-material fails. To generate such a field, in the material input file,
-a random quantity needs be added to the base value:
-\begin{cpp}
- sigma_c = $base$ uniform [$minimum$, $maximum$]
- sigma_c = $base$ weibull [$lambda$, $m$]
-\end{cpp}
-All parameters are real numbers. For the uniform distribution, minimum
-and maximum values have to be specified. The Weibull distribution is
-characterized by the following cumulative distribution function:
-\begin{equation}
- F(\sigma_\mathrm{c}) = 1- e^{-\left(
- \frac{\sigma_c-\sigma_\mathrm{c, min}}{\lambda} \right)^m}
-\end{equation}
-which depends on $\sigma_\mathrm{c, min}$, which is the minimum value,
-and $m$ and $\lambda$, which are the shape parameter and the scale
-parameter.
-
-The following sections describe the constitutive models implemented in
-\akantu. In Appendix~\ref{app:material-parameters} a summary of the
-parameters for all materials of \akantu is provided.
-\subsection{Elasticity}\index{Material!Elastic}
-The elastic law is a commonly used constitutive relationship that can
-be used for a wide range of engineering materials (\eg metals,
-concrete, rock, wood, glass, rubber, etc.) provided that the strains
-remain small (\ie small deformation and stress lower than yield
-strength). The linear elastic behavior is described by Hooke's law,
-which states that the stress is linearly proportional to the applied
-strain (material behaves like an ideal spring), as illustrated in
-Figure~\ref{fig:smm:cl:elastic}.
-\begin{figure}[!htb]
- \begin{center}
-
- \subfloat[]{
- \includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/stress_strain_el.pdf}
- \label{fig:smm:cl:elastic:stress_strain} }
- \hspace{0.05\textwidth} \subfloat[]{
- \raisebox{0.125\textwidth}{\includegraphics[width=0.25\textwidth,keepaspectratio=true]{figures/hooke_law.pdf}}
- \label{fig:smm:cl:elastic:hooke} }
- \caption{(a) Stress-strain curve for elastic material and (b)
- schematic representation of Hooke's law, denoted as a spring.}
- \label{fig:smm:cl:elastic}
- \end{center}
-\end{figure}
-The equation that relates the strains to the
-displacements is: % First the strain is computed (at every gauss
-point) from the displacements as follows:
-\begin{equation}
-\label{eqn:smm:strain_inf}
-\mat{\epsilon} =
-\frac{1}{2} \left[ \mat{F}-\mat{I}+(\mat{F}-\mat{I})^T \right]
-\end{equation}
-where $\mat{\epsilon}$ represents the infinitesimal strain tensor,
-$\mat{F} = \nabla_{\!\!\vec{X}}\vec{x}$ the deformation gradient
-tensor and $\mat{I}$ the identity matrix. The constitutive equation
-for isotropic homogeneous media can be expressed as:
-\begin{equation}
-\label{eqn:smm:material:constitutive_elastic}
-\mat{\sigma } =\lambda\mathrm{tr}(\mat{\epsilon})\mat{I}+2 \mu\mat{\epsilon}
-\end{equation}
-where $\mat{\sigma}$ is the Cauchy stress tensor
-($\lambda$ and $\mu$ are the the first and second Lame's
-coefficients). Besides the name and density, one has to specify the
-following properties in the material file: \code{E} (Young's modulus),
-\code{nu} (Poisson's ratio) and (for 2D) \code{Plane\_Stress} (if set
-to zero or not specified plane strain conditions are assumed for a
-plain analysis, otherwise plane stress conditions are applied).
-
-
\input{manual-constitutive-laws}
\section{Adding a New Constitutive Law}\index{Material!create a new
material}
There are several constitutive laws in \akantu as described in the
previous Section~\ref{sect:smm:CL}. It is also possible to use a
user-defined material for the simulation. These materials are referred
to as local materials since they are local to the example of the user
and not part of the \akantu library. To define a new local material,
two files (\code {material\_XXX.hh} and \code{material\_XXX.cc}) have
to be provided where \code{XXX} is the name of the new material. The
header file \code {material\_XXX.hh} defines the interface of your
custom material. Its implementation is provided in the
\code{material\_XXX.cc}. The new law must inherit from the
\code{Material} class or any other existing material class. It is
therefore necessary to include the interface of the parent material
in the header file of your local material and indicate the inheritance
in the declaration of the class:
\begin{cpp}
/* -------------------------------------------------------------------------- */
#include "material.hh"
/* -------------------------------------------------------------------------- */
#ifndef __AKANTU_MATERIAL_XXX_HH__
#define __AKANTU_MATERIAL_XXX_HH__
__BEGIN_AKANTU__
class MaterialXXX : public Material {
/// declare here the interface of your material
};
\end{cpp}
-In the header file the user also needs to declare all the members of
-the new material. These include the parameters that a read from the
-material input file, as well as any other material parameters that will be computed during the simulation and internal variables.\\
-In the following the example of a new damage material will be presented. In this case the parameters in the material will consist of the Young's modulus, the Poisson coefficient, the resistance to damage and the damage threshold. The material will then from these values compute its Lam\'{e} coefficients and its bulk modulus. Furthermore, the user has to add a new internal variable \code{damage} in order to store the amount of damage at each quadrature point in each step of the simulation. For this specific material the member declaration inside the class will look like follows:
+In the header file the user also needs to declare all the members of the new
+material. These include the parameters that a read from the
+material input file, as well as any other material parameters that will be computed during the simulation and internal variables.\\
+In the following the example of a new damage material will be presented. In this
+case the parameters in the material will consist of the Young's modulus, the
+Poisson coefficient, the resistance to damage and the damage threshold. The
+material will then from these values compute its Lam\'{e} coefficients and its
+bulk modulus. Furthermore, the user has to add a new internal variable
+\code{damage} in order to store the amount of damage at each quadrature point in
+each step of the simulation. For this specific material the member declaration
+inside the class will look like follows:
\begin{cpp}
class LocalMaterialDamage : public Material {
/// declare constructors/destructors here
/// declare methods and accessors here
/* ------------------------------------------------------------------------ */
/* Class Members */
/* ------------------------------------------------------------------------ */
AKANTU_GET_MACRO_BY_ELEMENT_TYPE_CONST(Damage, damage, Real);
private:
/// the young modulus
Real E;
/// Poisson coefficient
Real nu;
/// First Lamé coefficient
Real lambda;
/// Second Lamé coefficient (shear modulus)
Real mu;
/// resistance to damage
Real Yd;
/// damage threshold
Real Sd;
/// Bulk modulus
Real kpa;
/// damage internal variable
InternalField<Real> damage;
};
\end{cpp}
In order to enable to print the material parameters at any point in
the user's example file using the standard output stream by typing:
\begin{cpp}
for (UInt m = 0; m < model.getNbMaterials(); ++m)
std::cout << model.getMaterial(m) << std::endl;
\end{cpp}
the standard output stream operator has to be redefined. This should be done at the end of the header file:
\begin{cpp}
class LocalMaterialDamage : public Material {
/// declare here the interace of your material
}:
/* -------------------------------------------------------------------------- */
/* inline functions */
/* -------------------------------------------------------------------------- */
/// standard output stream operator
inline std::ostream & operator <<(std::ostream & stream, const LocalMaterialDamage & _this)
{
_this.printself(stream);
return stream;
}
\end{cpp}
However, the user still needs to register the material parameters that
should be printed out. The registration is done during the call of the
constructor. Like all definitions the implementation of the
constructor has to be written in the \code{material\_XXX.cc}
file. However, the declaration has to be provided in the
\code{material\_XXX.hh} file:
\begin{cpp}
class LocalMaterialDamage : public Material {
/* ------------------------------------------------------------------------ */
/* Constructors/Destructors */
/* ------------------------------------------------------------------------ */
public:
LocalMaterialDamage(SolidMechanicsModel & model, const ID & id = "");
};
\end{cpp}
The user can now define the implementation of the constructor in the
\code{material\_XXX.cc} file:
\begin{cpp}
/* -------------------------------------------------------------------------- */
#include "local_material_damage.hh"
#include "solid_mechanics_model.hh"
__BEGIN_AKANTU__
/* -------------------------------------------------------------------------- */
LocalMaterialDamage::LocalMaterialDamage(SolidMechanicsModel & model,
const ID & id) :
Material(model, id),
damage("damage", *this) {
AKANTU_DEBUG_IN();
this->registerParam("E" , E , 0. , _pat_parsable, "Young's modulus" );
this->registerParam("nu" , nu , 0.5 , _pat_parsable, "Poisson's ratio" );
this->registerParam("lambda" , lambda , _pat_readable, "First Lamé coefficient" );
this->registerParam("mu" , mu , _pat_readable, "Second Lamé coefficient");
this->registerParam("kapa" , kpa , _pat_readable, "Bulk coefficient" );
this->registerParam("Yd" , Yd , 50., _pat_parsmod);
this->registerParam("Sd" , Sd , 5000., _pat_parsmod);
damage.initialize(1);
AKANTU_DEBUG_OUT();
}
\end{cpp}
During the intializer list the reference to the model and the material id are assigned and the constructor of the internal field is called. Inside the scope of the constructor the internal values have to be initialized and the parameters, that should be printed out, are registered with the function:
\code{registerParam}\index{Material!registerParam}:
\begin{cpp}
void registerParam(name of the parameter (key in the material file),
member variable,
default value (optional parameter),
access permissions,
description);
\end{cpp}
The following table lists the all available access permissions:
\begin{center}
\begin{tabular}{ll}
\toprule
access permission & description\\
\midrule
\_pat\_internal & Parameter can only be output when the material is printed.\\
\_pat\_writable & User can write into the parameter. The parameter is output when the material is printed.\\
\_pat\_readable & User can read the parameter. \newline The parameter is output when the material is printed.\\
\_pat\_modifiable & Parameter is writable and readable.\\
\_pat\_parsable & Parameter can be parsed, \textit{i.e.} read from the input file.\\
\_pat\_parsmod & Parameter is modifiable and parsable.\\
\bottomrule
\end{tabular}
\end{center}
In order to implement the new constitutive law the user needs to
specify how the additional material parameters, that are not
defined in the input material file, should be calculated. Furthermore,
it has to be defined how stresses and the stable time step should be
computed for the new local material. In the case of implicit
simulations, in addition, the computation of the tangent stiffness needs
to be defined. Therefore, the user needs to redefine the following
functions of the parent material:
\begin{cpp}
void initMaterial();
// for explicit and implicit simulations void
computeStress(ElementType el_type, GhostType ghost_type = _not_ghost);
// for implicit simulations
void computeTangentStiffness(const ElementType & el_type,
Array<Real> & tangent_matrix,
GhostType ghost_type = _not_ghost);
// for explicit and implicit simulations
Real getStableTimeStep(Real h, const Element & element);
\end{cpp}
In the following a detailed description of these functions is provided:
\begin{itemize}
\item \code{initMaterial}:~ This method is called after the material
file is fully read and the elements corresponding to each material
are assigned. Some of the frequently used constant parameters are
calculated in this method. For example, the Lam\'{e} constants of
elastic materials can be considered as such parameters.
\item \code{computeStress}:~ In this method, the stresses are
computed based on the constitutive law as a function of the strains of the
quadrature points. For example, the stresses for the elastic
material are calculated based on the following formula:
\begin{equation}
\label{eqn:smm:constitutive_elastic}
\mat{\sigma } =\lambda\mathrm{tr}(\mat{\varepsilon})\mat{I}+2 \mu \mat{\varepsilon}
\end{equation}
Therefore, this method contains a loop on all quadrature points
assigned to the material using the
\code{MATERIAL\_STRESS\_QUADRATURE\_POINT\_LOOP\_BEGIN;} and
\code{MATERIAL\_STRESS\_QUADRATURE\_POINT\_LOOP\_END;} macros.
\begin{cpp}
MATERIAL_STRESS_QUADRATURE_POINT_LOOP_BEGIN;
// sigma <- f(grad_u)
MATERIAL_STRESS_QUADRATURE_POINT_LOOP_END;
\end{cpp}
\note{The strain vector in \akantu contains the values of $\nabla
\vec{u}$, i.e. it is really the \emph{displacement gradient},}
\item \code{computeTangentStiffness}:~ This method is called when
the tangent to the stress-strain curve is desired (see Fig \ref
{fig:smm:AL:K}). For example, it is called in the implicit solver
when the stiffness matrix for the regular elements is assembled
based on the following formula:
\begin{equation}
\label{eqn:smm:constitutive_elasc} \mat{K }
=\int{\mat{B^T}\mat{D(\varepsilon)}\mat{B}}
\end{equation}
Therefore, in this method, the \code{tangent} matrix (\mat{D}) is
computed for a given strain.
\note{ The \code{tangent} matrix is a $4^{th}$ order tensor which is
stored as a matrix in Voigt notation.}
\begin{figure}[!htb]
\begin{center}
\includegraphics[width=0.4\textwidth,keepaspectratio=true]{figures/tangent.pdf}
\caption{Tangent to the stress-strain curve.}
\label{fig:smm:AL:K}
\end{center}
\end{figure}
\item \code{getStableTimeStep}:~ The stable time step should be
calculated based on \eqref{eqn:smm:explicit:stabletime} for the
conditionally stable methods. This function depends on the longitudinal wave
speed which changes for different materials and strains. Therefore,
the value of this velocity should be defined in this method for
each new material.
\note {In case of need, the calculation of the longitudinal and
shear wave speeds can be done in separate functions
(\code{getPushWaveSpeed} and \code{getShearWaveSpeed}).
Therefore, the first function can be called for calculation of the
stable time step.}
\end{itemize}
Once the declaration and implementation of the new material has been
completed, this material can be used in the user's example by including the header file:
\begin{cpp}
#include "material_XXX.hh"
\end{cpp}
For existing materials, as mentioned in Section~\ref{sect:smm:CL}, by
default, the materials are initialized inside the method
\code{initFull}. If a local material should be used instead, the
initialization of the material has to be postponed until the local
material is registered in the model. Therefore, the model is
initialized with the boolean for skipping the material initialization
equal to true:
\begin{cpp}
/// model initialization
model.initFull(SolidMechanicsModelOptions(_explicit_lumped_mass,true));
\end{cpp}
Once the model has been initialized, the local material needs
to be registered in the model:
\begin{cpp}
model.registerNewCustomMaterials<XXX>("name_of_local_material");
\end{cpp}
Only at this point the material can be initialized:
\begin{cpp}
model.initMaterials();
\end{cpp}
A full example for adding a new damage law can be found in
\shellcode{\examplesdir/new\_material}.
%%% Local Variables: %%% mode: latex %%% TeX-master: "manual" %%% End:
diff --git a/doc/manual/manual.sty b/doc/manual/manual.sty
index 31f61d0c1..886b473d2 100644
--- a/doc/manual/manual.sty
+++ b/doc/manual/manual.sty
@@ -1,293 +1,309 @@
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% %%
%% LaTeX STYLE SHEET for AKANTU DOCUMENTATION %%
%% %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Geometry
\usepackage{a4wide}
\usepackage{geometry}
\geometry{
pdftex=true,
twoside=true,
margin=20mm,
bottom=20mm,
top=20mm,
bindingoffset=5.5mm
}
% Font encoding
\usepackage[T1]{fontenc}
\usepackage{palatino}
% Line spacing
\linespread{1.05}\selectfont
% Allow spaces to be added at the end of macro
\usepackage{xspace}
% Mathematics (including correct font encoding, after amsmath)
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{pxfonts}
\setlength\mathindent{2em}
% Some unit stuff
\usepackage[squaren]{SIunits}
% Easy tables
\usepackage{booktabs}
\usepackage{longtable}
\usepackage{rotating} % For sideways tables
\usepackage{multirow} % For larger cells
% A special environment for element tables
\newcommand\elemline{}%{\noalign{\smallskip}\hline\noalign{\smallskip}}
\newcommand\elemcooroned{\ensuremath{\left(\xi\right)}}
\newcommand\elemcoortwod{\ensuremath{\left(\xi\;,\;\eta\right)}}
\newcommand\elemcoorthreed{\ensuremath{\left(\xi\;,\;\eta\;,\;\zeta\right)}}
\newcommand\elemdshapeoned{\ensuremath{\left(\partial N_i/\partial\xi\right)}}
\newcommand\elemdshapetwod{\ensuremath{\left(\partial N_i/\partial\xi\;,\;\partial N_i/\partial\eta\right)}}
\newcommand\elemdshapethreed{\ensuremath{\left(\partial N_i/\partial\xi\;,\;\partial N_i/\partial\eta\;,\;\partial N_i/\partial\zeta\right)}}
\newcommand\inelemone[1]{\ensuremath{#1}}
%\newcommand\inelemtwo[2]{\ensuremath{\begin{pmatrix} {\; #1} & \!,\! & {#2 \;} \end{pmatrix}}}
%\newcommand\inelemthree[3]{\ensuremath{\begin{pmatrix} {\; #1} & \!,\! & {#2} & \!,\! & {#3 \;} \end{pmatrix}}}
\newcommand\inelemtwo[2]{\ensuremath{\left( {\; #1} \; , \; {#2 \;} \right)}}
\newcommand\inelemthree[3]{\ensuremath{\left( {\; #1} \; , \; {#2} \; , \; {#3 \;} \right)}}
\newcommand\inquadone[1]{\ensuremath{#1}}
\newcommand\inquadtwo[2]{\ensuremath{\left(\, #1 \, , \, #2 \,\right)}}
\newcommand\inquadthree[3]{\ensuremath{\left(\, #1 \, , \, #2 \, , \, #3 \,\right)}}
%\newcommand\quada{\tfrac{1}{20}\left(5-\sqrt{5}\right)}
%\newcommand\quadb{\tfrac{1}{20}\left(5+3\sqrt{5}\right)}
\newcommand\quada{\tfrac{\left(5-\sqrt{5}\right)}{20}}
\newcommand\quadb{\tfrac{\left(5+3\sqrt{5}\right)}{20}}
\newenvironment{Element}[1]
{%\begin{table*}[!htbp]
\footnotesize
\ifthenelse{\equal{#1}{1D}}{\renewcommand{\arraystretch}{1.50}}{}
\ifthenelse{\equal{#1}{2D}}{\renewcommand{\arraystretch}{1.60}}{}
\ifthenelse{\equal{#1}{3D}}{\renewcommand{\arraystretch}{1.70}}{}
\textbf{{\normalsize Element properties}}\\%\vspace*{0.5\baselineskip}
\begin{tabular}{cccc}\\[-2ex]
\toprule
% \multicolumn{4}{l}{\textbf{{\normalsize Element properties}}} \\
Node ($i$) & Coord.
\ifthenelse{\equal{#1}{1D}}{\elemcooroned}{}
\ifthenelse{\equal{#1}{2D}}{\elemcoortwod}{}
\ifthenelse{\equal{#1}{3D}}{\elemcoorthreed}{}
& Shape function ($N_{i}$)
& {Derivative
\ifthenelse{\equal{#1}{1D}}{\elemdshapeoned}{}
\ifthenelse{\equal{#1}{2D}}{\elemdshapetwod}{}
\ifthenelse{\equal{#1}{3D}}{\elemdshapethreed}{}
}\\
\midrule
\elemline
}
{\bottomrule
\end{tabular}%\end{table*}
}
\newenvironment{QuadPoints}[1]
{\vspace*{\baselineskip}
%\begin{table*}[!htbp]
\footnotesize
\renewcommand{\arraystretch}{1.50}
\noindent\textbf{{\normalsize Gaussian quadrature points}}\newline%\vspace*{0.5\baselineskip}
\begin{tabular}{#1}\\[-2ex]
\toprule
}
{\bottomrule\end{tabular}}%\end{table*}}
+\usepackage{xparse}
+
+\NewDocumentEnvironment{MaterialDesc}{m m}{
+ Keyword: \code{#1}\\[.8ex]
+ \noindent Description here: \ref{#2}\\[0.8ex]
+ \noindent Parameters:\\[-4ex]
+ \begin{itemize}
+ \setlength{\topsep}{0ex}%
+ \setlength{\parsep}{0cm}%
+ \setlength{\itemsep}{0cm}%
+ \setlength{\parskip}{0cm}%
+ }{\end{itemize}}
+
+\newcommand{\matparam}[3]{\item \code{#1}: (\emph{#2}) #3}
+\newcommand{\matinherit}[1]{\usebox{#1}}
+
% Nice coloring
\usepackage[dvipsnames,usenames,table]{xcolor}
\definecolor{RED}{rgb}{1,0,0}
\definecolor{cppbg}{HTML}{EBF2F2}
\definecolor{shellbg}{HTML}{F5EDE4}
\definecolor{commentcolor}{HTML}{101280}
% Allow for the use of listings
\usepackage{listings}
% Create an index
\usepackage{makeidx}
% Figure handling
\usepackage{graphics}
\usepackage{epsfig}
\usepackage[lofdepth,lotdepth]{subfig}
\usepackage{tikz}
\usetikzlibrary{decorations}
\renewcommand{\floatpagefraction}{.6} % default: .5
\renewcommand\topfraction{0.9} % 90% of page top can be a float (Standard 0.7)
\renewcommand\bottomfraction{0.1} % 10% of page bottom can be a float (Standard 0.3)
\renewcommand\textfraction{0.1} % only 10% of page must to be text (Standard 0.2)
% Removes parenthese around subfig number
\renewcommand*{\thesubfigure}{\alph{subfigure}}
% Create a new list style for C++
\lstdefinestyle{C++}{
language=C++, % the language of the code
basicstyle=\small\ttfamily, % Without beramono, we'd get cmtt, the teletype font.
commentstyle=\color{commentcolor}\itshape,
keywordstyle=\color{DarkOrchid}\bfseries,
% fontadjust,
% numbers=left, % where to put the line-numbers
% numberstyle=\tiny, % the size of the fonts that are used for the line-numbers
% stepnumber=2, % the step between two line-numbers. If it's 1, each line will
% be numbered
% numbersep=5pt, % how far the line-numbers are from the code
% showspaces=false, % show spaces adding particular underscores
showstringspaces=false, % underline spaces within strings
% showtabs=false, % show tabs within strings adding particular underscores
% frame=llines, % adds a frame around the code
% frame=tb,
tabsize=2, % sets default tabsize to 2 spaces
captionpos=b, % sets the caption-position to bottom
breaklines=true, % sets automatic line breaking
breakatwhitespace=false, % sets if automatic breaks should only happen at
% whitespace
% title=\lstname, % show the filename of files included with \lstinputlisting;
% also try caption instead of title
% escapeinside={\%*}{*)}, % if you want to add a comment within your code
xleftmargin=1cm,
xrightmargin=1cm,
mathescape=true,
escapechar=\%,
morekeywords={Real, UInt, Int},
columns=flexible,
keepspaces=true,
backgroundcolor=\color{cppbg}
}
% Create new list style for the shell
\lstdefinestyle{shell}{
language=bash, % the language of the code
basicstyle=\scriptsize\ttfamily, % Without beramono, we'd get cmtt, the teletype font.
showstringspaces=false, % underline spaces within strings
tabsize=2, % sets default tabsize to 2 spaces
captionpos=b, % sets the caption-position to bottom
breaklines=true, % sets automatic line breaking
breakatwhitespace=false,
xleftmargin=1cm,
xrightmargin=1cm,
escapechar=\%,
morekeywords={mkdir, make, ccmake, cmake},
columns=flexible,
keepspaces=true,
backgroundcolor=\color{shellbg}
}
% Set some derived listing environments
\lstnewenvironment{cpp}{\lstset{style=C++}}{}
\lstnewenvironment{command}{\lstset{style=shell}}{}
% Make sure outputspace is white
\makeatletter
\def\lst@outputspace{{\ifx\lst@bkgcolor\empty\color{white}\else\lst@bkgcolor\fi\lst@visiblespace}}
\makeatother
% Renow a label in the itemized lists
\renewcommand{\labelitemi}{$\mathbf{\circ}$}
% Don't care so much about overfull h-boxes
\sloppy
% Penalty adjusments
%\widowpenalty=10000 % Single lines/word on beginning of page
%\clubpenalty=10000 % Single lines/word at end of page
%\hyphenpenalty=2000 % Hyphenate words
%\tolerance=250 % To adjust the hyphenation of words, increase the tolerance to discourage hyphenation
% the higher the value, the more ugly the gaps between words
% default \tolerance=200
%\hfuzz=10000pt % threshold when an overfull hbox is reported, default \hfuzz=0.1pt
%\vfuzz=10000pt % threshold to report an overfull vbox, default \vfuzz=0.1pt
%\hbadness=10000 % threshold to report an underfull \hbox
%\vbadness=10000 % threshold to report an underfull \vbox protokolliert wird.
\emergencystretch=0pt % causes a third attempt to fix bad paragraphs and defines a maximum limit to stretch them
% Insert an empty or a blank page
\newcommand{\insertemptypage}{\newpage\hbox{}\newpage}
\newcommand{\insertblankpage}{\newpage\thispagestyle{empty}\hbox{}\newpage}
% No page number on an empty page
\let\origdoublepage\cleardoublepage
\newcommand{\clearemptydoublepage}{%
\clearpage
{\thispagestyle{empty}\origdoublepage}%
}
\let\cleardoublepage\clearemptydoublepage
% A new ruler for in chapters
\newcommand\InChapterRule{\addvspace{\baselineskip}\rule{0.3\linewidth}{0.25pt}}
% New footnote style
%\def\@fnsymbol#1{\ifcase#1\or *\or \dagger\or \ddagger\or \mathchar "278\or \mathchar "27B\or \|\or **\or \dagger\dagger \or \ddagger\ddagger \else\@ctrerr\fi\relax}
\def\@fnsymbol#1{*\xspace\relax}
\renewcommand{\thefootnote}{\fnsymbol{footnote}} % Symbols rather than numbers
\def\footnoterule{\vspace*{0.5\baselineskip}\InChapterRule\vspace*{0.25\baselineskip}}
% Improved look of the Table of Contents
\usepackage[nottoc,notbib]{tocbibind}
\usepackage[dotinlabels]{titletoc}
\titlecontents{chapter}[1.4pc]
{\addvspace{0.6pc}\large\bfseries\filright}
{\contentslabel[\thecontentslabel.]{1.4pc}}
{\hspace{-1.4pc}}
{\hfill\contentspage}
[\addvspace{2pt}]
\titlecontents{section}[3.4pc]
{\filright}
{\contentslabel[\thecontentslabel]{2pc}}
{\hspace{-2pc}}
{\titlerule*[6pt]{.}\contentspage}
[]
\titlecontents{subsection}[5.0pc]
{\filright}
{\contentslabel[\thecontentslabel]{2.4pc}}
{}
{\titlerule*[6pt]{.}\contentspage}
[]
\setcounter{tocdepth}{2}
\newcommand\addspaceintoc{\addtocontents{toc}{\protect\addvspace{20pt}}}
% Change the appearance of the bibliography
\bibliographystyle{manual-bibliographystyle}
\renewcommand\bibname{References}
\usepackage{cite} % To sort citations: [2,10-14]
\let\oldthebibliography=\thebibliography
\let\endoldthebibliography=\endthebibliography
\renewenvironment{thebibliography}[1]
{ \begin{oldthebibliography}{#1}
% \small
\addcontentsline{toc}{chapter}{\bibname}
\setlength{\labelsep}{2mm}
\setlength{\parskip}{0\baselineskip}
\setlength{\itemsep}{0.24\baselineskip}
}
{ \end{oldthebibliography} }
% Hyperref
\usepackage{url}
\usepackage[pdftex,
bookmarks=true,
bookmarksnumbered=true,
% linkbordercolor={1 1 1},
% pdfborder={0 0 0},
pdfpagemode=UseOutlines
]{hyperref}
\hypersetup{
pdfauthor={Computational Solid Mechanics Laboratory - EPFL},
pdftitle={Akantu User's Guide},
pdfsubject={Open Source Finite Element Code - Akantu}
}
\newcommand\akantuPackageSection[2]{\paragraph*{\label{#2}#1}}
%\newcommand\akantuSubpackageName[3]{\item \index{Packages!#1}\hyperref[#2]{#1}}
%\newcommand\akantuPackageDependencyHead{~\\The dependencies are: \begin{itemize}}
%\newcommand\akantuPackageDependencyTail{\end{itemize}}
\newcommand\akantuSubpackageName[3]{ \textit{\index{Packages!#1}\hyperref[#2]{#1}}}
\newcommand\akantuPackageDependencyHead{~\\The dependencies are: }
\newcommand\akantuPackageDependencyTail{~\\}
diff --git a/packages/00_core.cmake b/packages/00_core.cmake
index af097463e..b328f43a8 100644
--- a/packages/00_core.cmake
+++ b/packages/00_core.cmake
@@ -1,405 +1,405 @@
#===============================================================================
# @file core.cmake
#
# @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Mon Nov 21 18:19:15 2011
#
# @brief package description for core
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
set(AKANTU_CORE ON CACHE INTERNAL "core package for Akantu" FORCE)
set(AKANTU_CORE_FILES
common/aka_array.cc
common/aka_array.hh
common/aka_array_tmpl.hh
common/aka_blas_lapack.hh
common/aka_circular_array.hh
common/aka_circular_array_inline_impl.cc
common/aka_common.cc
common/aka_common.hh
common/aka_common_inline_impl.cc
common/aka_csr.hh
common/aka_error.cc
common/aka_error.hh
common/aka_event_handler_manager.hh
common/aka_extern.cc
common/aka_fwd.hh
common/aka_grid_dynamic.hh
common/aka_math.cc
common/aka_math.hh
common/aka_math_tmpl.hh
common/aka_memory.cc
common/aka_memory.hh
common/aka_memory_inline_impl.cc
common/aka_random_generator.hh
common/aka_safe_enum.hh
common/aka_static_memory.cc
common/aka_static_memory.hh
common/aka_static_memory_inline_impl.cc
common/aka_static_memory_tmpl.hh
common/aka_typelist.hh
common/aka_types.hh
common/aka_visitor.hh
common/aka_voigthelper.hh
common/aka_voigthelper.cc
- fem/element_class.cc
- fem/element_class.hh
- fem/element_class_tmpl.hh
- fem/element_classes/element_class_hexahedron_8_inline_impl.cc
- fem/element_classes/element_class_pentahedron_6_inline_impl.cc
- fem/element_classes/element_class_point_1_inline_impl.cc
- fem/element_classes/element_class_quadrangle_4_inline_impl.cc
- fem/element_classes/element_class_quadrangle_8_inline_impl.cc
- fem/element_classes/element_class_segment_2_inline_impl.cc
- fem/element_classes/element_class_segment_3_inline_impl.cc
- fem/element_classes/element_class_tetrahedron_10_inline_impl.cc
- fem/element_classes/element_class_tetrahedron_4_inline_impl.cc
- fem/element_classes/element_class_triangle_3_inline_impl.cc
- fem/element_classes/element_class_triangle_6_inline_impl.cc
+ fe_engine/element_class.cc
+ fe_engine/element_class.hh
+ fe_engine/element_class_tmpl.hh
+ fe_engine/element_classes/element_class_hexahedron_8_inline_impl.cc
+ fe_engine/element_classes/element_class_pentahedron_6_inline_impl.cc
+ fe_engine/element_classes/element_class_point_1_inline_impl.cc
+ fe_engine/element_classes/element_class_quadrangle_4_inline_impl.cc
+ fe_engine/element_classes/element_class_quadrangle_8_inline_impl.cc
+ fe_engine/element_classes/element_class_segment_2_inline_impl.cc
+ fe_engine/element_classes/element_class_segment_3_inline_impl.cc
+ fe_engine/element_classes/element_class_tetrahedron_10_inline_impl.cc
+ fe_engine/element_classes/element_class_tetrahedron_4_inline_impl.cc
+ fe_engine/element_classes/element_class_triangle_3_inline_impl.cc
+ fe_engine/element_classes/element_class_triangle_6_inline_impl.cc
- fem/fe_engine.cc
- fem/fe_engine.hh
- fem/fe_engine_inline_impl.cc
- fem/fe_engine_template.hh
- fem/fe_engine_template_tmpl.hh
- fem/geometrical_data_tmpl.hh
- fem/geometrical_element.cc
- fem/integration_element.cc
- fem/integrator.hh
- fem/integrator_gauss.hh
- fem/integrator_gauss_inline_impl.cc
- fem/interpolation_element.cc
- fem/interpolation_element_tmpl.hh
- fem/shape_functions.hh
- fem/shape_functions_inline_impl.cc
- fem/shape_lagrange.cc
- fem/shape_lagrange.hh
- fem/shape_lagrange_inline_impl.cc
- fem/shape_linked.cc
- fem/shape_linked.hh
- fem/shape_linked_inline_impl.cc
+ fe_engine/fe_engine.cc
+ fe_engine/fe_engine.hh
+ fe_engine/fe_engine_inline_impl.cc
+ fe_engine/fe_engine_template.hh
+ fe_engine/fe_engine_template_tmpl.hh
+ fe_engine/geometrical_data_tmpl.hh
+ fe_engine/geometrical_element.cc
+ fe_engine/integration_element.cc
+ fe_engine/integrator.hh
+ fe_engine/integrator_gauss.hh
+ fe_engine/integrator_gauss_inline_impl.cc
+ fe_engine/interpolation_element.cc
+ fe_engine/interpolation_element_tmpl.hh
+ fe_engine/shape_functions.hh
+ fe_engine/shape_functions_inline_impl.cc
+ fe_engine/shape_lagrange.cc
+ fe_engine/shape_lagrange.hh
+ fe_engine/shape_lagrange_inline_impl.cc
+ fe_engine/shape_linked.cc
+ fe_engine/shape_linked.hh
+ fe_engine/shape_linked_inline_impl.cc
io/dumper/dumpable.hh
io/dumper/dumpable_inline_impl.hh
io/mesh_io.cc
io/mesh_io.hh
io/mesh_io/mesh_io_abaqus.cc
io/mesh_io/mesh_io_abaqus.hh
io/mesh_io/mesh_io_diana.cc
io/mesh_io/mesh_io_diana.hh
io/mesh_io/mesh_io_msh.cc
io/mesh_io/mesh_io_msh.hh
io/model_io.cc
io/model_io.hh
io/parser/algebraic_parser.hh
io/parser/input_file_parser.hh
io/parser/parsable.cc
io/parser/parsable.hh
io/parser/parsable_tmpl.hh
io/parser/parser.cc
io/parser/parser.hh
io/parser/parser_tmpl.hh
io/parser/cppargparse/cppargparse.hh
io/parser/cppargparse/cppargparse.cc
io/parser/cppargparse/cppargparse_tmpl.hh
mesh/element_group.cc
mesh/element_group.hh
mesh/element_group_inline_impl.cc
mesh/element_type_map.hh
mesh/element_type_map_tmpl.hh
mesh/group_manager.cc
mesh/group_manager.hh
mesh/group_manager_inline_impl.cc
mesh/mesh.cc
mesh/mesh.hh
mesh/mesh_filter.hh
mesh/mesh_data.cc
mesh/mesh_data.hh
mesh/mesh_data_tmpl.hh
mesh/mesh_inline_impl.cc
mesh/node_group.cc
mesh/node_group.hh
mesh/node_group_inline_impl.cc
mesh_utils/mesh_partition.cc
mesh_utils/mesh_partition.hh
mesh_utils/mesh_partition/mesh_partition_mesh_data.cc
mesh_utils/mesh_partition/mesh_partition_mesh_data.hh
mesh_utils/mesh_partition/mesh_partition_scotch.hh
mesh_utils/mesh_pbc.cc
mesh_utils/mesh_utils.cc
mesh_utils/mesh_utils.hh
mesh_utils/mesh_utils_inline_impl.cc
mesh_utils/mesh_graph.cc
mesh_utils/mesh_graph.hh
model/boundary_condition.hh
model/boundary_condition_functor.hh
model/boundary_condition_functor_inline_impl.cc
model/boundary_condition_tmpl.hh
model/integration_scheme/generalized_trapezoidal.hh
model/integration_scheme/generalized_trapezoidal_inline_impl.cc
model/integration_scheme/integration_scheme_1st_order.hh
model/integration_scheme/integration_scheme_2nd_order.hh
model/integration_scheme/newmark-beta.hh
model/integration_scheme/newmark-beta_inline_impl.cc
model/model.cc
model/model.hh
model/model_inline_impl.cc
model/solid_mechanics/material.cc
model/solid_mechanics/material.hh
model/solid_mechanics/material_inline_impl.cc
model/solid_mechanics/material_list.hh
model/solid_mechanics/material_random_internal.hh
model/solid_mechanics/material_selector.hh
model/solid_mechanics/material_selector_tmpl.hh
model/solid_mechanics/materials/internal_field.hh
model/solid_mechanics/materials/internal_field_tmpl.hh
model/solid_mechanics/materials/material_elastic.cc
model/solid_mechanics/materials/material_elastic.hh
model/solid_mechanics/materials/material_elastic_inline_impl.cc
model/solid_mechanics/materials/material_thermal.cc
model/solid_mechanics/materials/material_thermal.hh
model/solid_mechanics/materials/random_internal_field.hh
model/solid_mechanics/materials/random_internal_field_tmpl.hh
model/solid_mechanics/solid_mechanics_model.cc
model/solid_mechanics/solid_mechanics_model.hh
model/solid_mechanics/solid_mechanics_model_inline_impl.cc
model/solid_mechanics/solid_mechanics_model_mass.cc
model/solid_mechanics/solid_mechanics_model_material.cc
model/solid_mechanics/solid_mechanics_model_tmpl.hh
model/solid_mechanics/solid_mechanics_model_event_handler.hh
model/solid_mechanics/materials/material_damage/material_damage.hh
model/solid_mechanics/materials/material_damage/material_damage_tmpl.hh
model/solid_mechanics/materials/material_damage/material_marigo.cc
model/solid_mechanics/materials/material_damage/material_marigo.hh
model/solid_mechanics/materials/material_damage/material_marigo_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_mazars.cc
model/solid_mechanics/materials/material_damage/material_mazars.hh
model/solid_mechanics/materials/material_damage/material_mazars_inline_impl.cc
model/solid_mechanics/materials/material_finite_deformation/material_neohookean.cc
model/solid_mechanics/materials/material_finite_deformation/material_neohookean.hh
model/solid_mechanics/materials/material_finite_deformation/material_neohookean_inline_impl.cc
model/solid_mechanics/materials/material_plastic/material_plastic.cc
model/solid_mechanics/materials/material_plastic/material_plastic.hh
model/solid_mechanics/materials/material_plastic/material_plastic_inline_impl.cc
model/solid_mechanics/materials/material_plastic/material_linear_isotropic_hardening.cc
model/solid_mechanics/materials/material_plastic/material_linear_isotropic_hardening.hh
model/solid_mechanics/materials/material_plastic/material_linear_isotropic_hardening_inline_impl.cc
model/solid_mechanics/materials/material_viscoelastic/material_standard_linear_solid_deviatoric.cc
model/solid_mechanics/materials/material_viscoelastic/material_standard_linear_solid_deviatoric.hh
solver/solver.cc
solver/solver.hh
solver/sparse_matrix.cc
solver/sparse_matrix.hh
solver/sparse_matrix_inline_impl.cc
synchronizer/communication_buffer.hh
synchronizer/communication_buffer_inline_impl.cc
synchronizer/data_accessor.cc
synchronizer/data_accessor.hh
synchronizer/data_accessor_inline_impl.cc
synchronizer/distributed_synchronizer.cc
synchronizer/distributed_synchronizer.hh
synchronizer/distributed_synchronizer_tmpl.hh
synchronizer/dof_synchronizer.cc
synchronizer/dof_synchronizer.hh
synchronizer/dof_synchronizer_inline_impl.cc
synchronizer/filtered_synchronizer.cc
synchronizer/filtered_synchronizer.hh
synchronizer/mpi_type_wrapper.hh
synchronizer/pbc_synchronizer.cc
synchronizer/pbc_synchronizer.hh
synchronizer/real_static_communicator.hh
synchronizer/static_communicator.cc
synchronizer/static_communicator.hh
synchronizer/static_communicator_dummy.hh
synchronizer/static_communicator_inline_impl.hh
synchronizer/synchronizer.cc
synchronizer/synchronizer.hh
synchronizer/synchronizer_registry.cc
synchronizer/synchronizer_registry.hh
)
set(AKANTU_CORE_DEB_DEPEND
libboost-dev
)
set(AKANTU_CORE_TESTS
test_csr
test_facet_element_mapping
test_facet_extraction_tetrahedron_4
test_facet_extraction_triangle_3
test_grid
test_interpolate_stress
test_local_material
test_material_damage_non_local
test_material_thermal
test_matrix
test_mesh_boundary
test_mesh_data
test_mesh_io_msh
test_mesh_io_msh_physical_names
test_mesh_partitionate_mesh_data
test_parser
test_pbc_tweak
test_purify_mesh
test_solid_mechanics_model_bar_traction2d
test_solid_mechanics_model_bar_traction2d_structured
test_solid_mechanics_model_bar_traction2d_structured_pbc
test_solid_mechanics_model_boundary_condition
test_solid_mechanics_model_circle_2
test_solid_mechanics_model_cube3d
test_solid_mechanics_model_cube3d_pbc
test_solid_mechanics_model_cube3d_tetra10
test_solid_mechanics_model_square
test_solid_mechanics_model_material_eigenstrain
test_static_memory
test_surface_extraction_tetrahedron_4
test_surface_extraction_triangle_3
test_vector
test_vector_iterator
test_weight
test_math
test_material_standard_linear_solid_deviatoric_relaxation
test_material_standard_linear_solid_deviatoric_relaxation_tension
test_material_plasticity
)
set(AKANTU_CORE_MANUAL_FILES
manual.sty
manual.cls
manual.tex
manual-macros.sty
manual-titlepages.tex
manual-introduction.tex
manual-gettingstarted.tex
manual-io.tex
manual-solidmechanicsmodel.tex
manual-constitutive-laws.tex
manual-lumping.tex
manual-elements.tex
manual-appendix-elements.tex
manual-appendix-materials.tex
manual-backmatter.tex
manual-bibliography.bib
manual-bibliographystyle.bst
figures/bc_and_ic_example.pdf
figures/boundary.pdf
figures/boundary.svg
figures/dirichlet.pdf
figures/dirichlet.svg
figures/doc_wheel.pdf
figures/doc_wheel.svg
figures/dynamic_analysis.png
figures/explicit_dynamic.pdf
figures/explicit_dynamic.svg
figures/static.pdf
figures/static.svg
figures/hooke_law.pdf
figures/hot-point-1.png
figures/hot-point-2.png
figures/implicit_dynamic.pdf
figures/implicit_dynamic.svg
figures/insertion.pdf
figures/interpolate.pdf
figures/interpolate.svg
figures/problemDomain.pdf_tex
figures/problemDomain.pdf
figures/static_analysis.png
figures/stress_strain_el.pdf
figures/tangent.pdf
figures/tangent.svg
figures/vectors.pdf
figures/vectors.svg
figures/stress_strain_neo.pdf
figures/visco_elastic_law.pdf
figures/isotropic_hardening_plasticity.pdf
figures/elements/hexahedron_8.pdf
figures/elements/hexahedron_8.svg
figures/elements/quadrangle_4.pdf
figures/elements/quadrangle_4.svg
figures/elements/quadrangle_8.pdf
figures/elements/quadrangle_8.svg
figures/elements/segment_2.pdf
figures/elements/segment_2.svg
figures/elements/segment_3.pdf
figures/elements/segment_3.svg
figures/elements/tetrahedron_10.pdf
figures/elements/tetrahedron_10.svg
figures/elements/tetrahedron_4.pdf
figures/elements/tetrahedron_4.svg
figures/elements/triangle_3.pdf
figures/elements/triangle_3.svg
figures/elements/triangle_6.pdf
figures/elements/triangle_6.svg
figures/elements/xtemp.pdf
)
find_program(READLINK_COMMAND readlink)
find_program(ADDR2LINE_COMMAND addr2line)
mark_as_advanced(READLINK_COMMAND)
mark_as_advanced(ADDR2LINE_COMMAND)
include(CheckFunctionExists)
check_function_exists(clock_gettime _clock_gettime)
if(NOT _clock_gettime)
set(AKANTU_USE_OBSOLETE_GETTIMEOFDAY ON)
else()
set(AKANTU_USE_OBSOLETE_GETTIMEOFDAY OFF)
endif()
list(APPEND AKANTU_BOOST_COMPONENTS
graph
)
set(AKANTU_CORE_DOCUMENTATION
"
This package is the core engine of Akantu. It depends on:
\\begin{itemize}
\\item A C++ compiler (\\href{http://gcc.gnu.org/}{GCC} >= 4, or \\href{https://software.intel.com/en-us/intel-compilers}{Intel}).
\\item The cross-platform, open-source \\href{http://www.cmake.org/}{CMake} build system.
\\item The \\href{http://www.boost.org/}{Boost} C++ portable libraries.
\\item The \\href{http://www.zlib.net/}{zlib} compression library.
\\end{itemize}
Under Ubuntu (14.04 LTS) the installation can be performed using the commands:
\\begin{command}
> sudo apt-get install build-essential cmake-curses-gui libboost-dev zlib1g-dev
\\end{command}
")
\ No newline at end of file
diff --git a/packages/05_igfem.cmake b/packages/05_igfem.cmake
index 3db136f18..e07545e58 100644
--- a/packages/05_igfem.cmake
+++ b/packages/05_igfem.cmake
@@ -1,47 +1,46 @@
#===============================================================================
# @file igfem.cmake
#
# @author Aurelia Cuba Ramos <aurelia.cubaramos@epfl.ch>
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Fri Mai 23 18:19:15 2011
#
# @brief package description for interface-enriched generalized IGFEM
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
option(AKANTU_IGFEM "Use Interface-enriched generalized FEM" OFF)
set(AKANTU_IGFEM_FILES
- fem/element_class_igfem.hh
- fem/shape_igfem.hh
- fem/shape_igfem_inline_impl.cc
- fem/element_class_igfem.cc
- fem/element_classes/element_class_igfem_triangle_3_inline_impl.cc
- fem/igfem_element.hh
- fem/igfem_element.cc
-
+ fe_engine/element_class_igfem.hh
+ fe_engine/shape_igfem.hh
+ fe_engine/shape_igfem_inline_impl.cc
+ fe_engine/element_class_igfem.cc
+ fe_engine/element_classes/element_class_igfem_triangle_3_inline_impl.cc
+ fe_engine/igfem_element.hh
+ fe_engine/igfem_element.cc
)
set(AKANTU_IGFEM_TESTS
test_igfem_integrate
)
diff --git a/packages/10_structural_mechanics.cmake b/packages/10_structural_mechanics.cmake
index 1515e8e11..e0dab4bee 100644
--- a/packages/10_structural_mechanics.cmake
+++ b/packages/10_structural_mechanics.cmake
@@ -1,74 +1,74 @@
#===============================================================================
# @file structural_mechanics.cmake
#
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Mon Nov 21 18:19:15 2011
#
# @brief package description for structural mechanics
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
option(AKANTU_STRUCTURAL_MECHANICS "Use Structural mechanics model package of Akantu" OFF)
add_internal_package_dependencies(STRUCTURAL_MECHANICS IMPLICIT)
set(AKANTU_STRUCTURAL_MECHANICS_FILES
- fem/element_class_structural.hh
- fem/element_classes/element_class_bernoulli_beam_inline_impl.cc
+ fe_engine/element_class_structural.hh
+ fe_engine/element_classes/element_class_bernoulli_beam_inline_impl.cc
io/mesh_io/mesh_io_msh_struct.cc
io/mesh_io/mesh_io_msh_struct.hh
io/model_io/model_io_ibarras.cc
io/model_io/model_io_ibarras.hh
model/structural_mechanics/structural_mechanics_model.cc
model/structural_mechanics/structural_mechanics_model.hh
model/structural_mechanics/structural_mechanics_model_boundary.cc
model/structural_mechanics/structural_mechanics_model_inline_impl.cc
)
set(AKANTU_STRUCTURAL_MECHANICS_DOC
manual/manual-structuralmechanicsmodel.tex
)
set(AKANTU_STRUCTURAL_MECHANICS_TESTS
test_structural_mechanics_model_bernoulli_beam_2
test_structural_mechanics_model_boundary_bernoulli_beam_2
test_structural_mechanics_model_bernoulli_beam_2_exemple_1_1
test_structural_mechanics_model_bernoulli_beam_2_complicated
test_structural_mechanics_model_bernoulli_beam_2_exemple_1_1_y
test_structural_mechanics_model_bernoulli_beam_3_exemple_1_1_xy
test_structural_mechanics_model_bernoulli_beam_3_exemple_1_1_zy
test_structural_mechanics_model_bernoulli_beam_3_local_force
test_structural_mechanics_model_bernoulli_beam_3_exercice_12_10_13
)
set(AKANTU_STRUCTURAL_MECHANICS_MANUAL_FILES
manual-structuralmechanicsmodel.tex
manual-structuralmechanicsmodel-elements.tex
figures/beam_example.pdf
figures/elements/bernoulli_2.pdf
figures/elements/bernoulli_2.svg
)
set(AKANTU_STRUCTURAL_MECHANICS_DOCUMENTATION "
This package activates the compilation for the Structural Mechanics engine of Akantu
")
diff --git a/packages/20_cohesive_element.cmake b/packages/20_cohesive_element.cmake
index 853c98a61..237ffae7b 100644
--- a/packages/20_cohesive_element.cmake
+++ b/packages/20_cohesive_element.cmake
@@ -1,111 +1,112 @@
#===============================================================================
# @file cohesive_element.cmake
#
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Tue Oct 16 14:05:02 2012
#
# @brief package description for cohesive elements
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
option(AKANTU_COHESIVE_ELEMENT "Use cohesive_element package of Akantu" OFF)
add_external_package_dependencies(cohesive_element lapack)
set(AKANTU_COHESIVE_ELEMENT_FILES
model/solid_mechanics/materials/material_cohesive_includes.hh
mesh_utils/cohesive_element_inserter.hh
mesh_utils/cohesive_element_inserter.cc
- fem/cohesive_element.cc
- fem/shape_cohesive.hh
- fem/cohesive_element.hh
- fem/fe_engine_template_cohesive.cc
+ fe_engine/cohesive_element.cc
+ fe_engine/shape_cohesive.hh
+ fe_engine/cohesive_element.hh
+ fe_engine/fe_engine_template_cohesive.cc
- fem/shape_cohesive_inline_impl.cc
+ fe_engine/shape_cohesive_inline_impl.cc
model/solid_mechanics/materials/material_cohesive/cohesive_internal_field_tmpl.hh
model/solid_mechanics/materials/material_cohesive/cohesive_internal_field.hh
model/solid_mechanics/materials/material_cohesive/material_cohesive_inline_impl.cc
model/solid_mechanics/solid_mechanics_model_cohesive.cc
model/solid_mechanics/fragment_manager.cc
model/solid_mechanics/materials/material_cohesive/material_cohesive.cc
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_linear.cc
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_bilinear.cc
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_exponential.cc
model/solid_mechanics/solid_mechanics_model_cohesive.hh
model/solid_mechanics/fragment_manager.hh
model/solid_mechanics/materials/material_cohesive/material_cohesive.hh
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_bilinear.hh
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_linear.hh
model/solid_mechanics/materials/material_cohesive/constitutive_laws/material_cohesive_exponential.hh
io/dumper/dumper_iohelper_tmpl_connectivity_field.hh
)
set(AKANTU_COHESIVE_ELEMENT_TESTS
test_cohesive_buildfacets_tetrahedron
test_cohesive_buildfacets_hexahedron
test_cohesive_intrinsic
test_cohesive_intrinsic_quadrangle
test_cohesive_extrinsic
test_cohesive_extrinsic_quadrangle
test_cohesive_buildfragments
test_cohesive_intrinsic_impl
)
set(AKANTU_COHESIVE_ELEMENT_DOC
manual/manual-cohesive_elements.tex
manual/manual-cohesive_elements_insertion.tex
manual/manual-cohesive_laws.tex
+ manual-appendix-materials-cohesive.tex
)
set(AKANTU_COHESIVE_ELEMENT_MANUAL_FILES
manual-cohesive_elements.tex
manual-cohesive_elements_insertion.tex
manual-cohesive_laws.tex
figures/cohesive2d.pdf
figures/cohesive3d.pdf
figures/cohesive_exponential.pdf
figures/linear_cohesive_law.pdf
figures/bilinear_cohesive_law.pdf
)
set(AKANTU_COHESIVE_ELEMENT_DOCUMENTATION
"This package activates the cohesive elements engine within Akantu.
It depends on:
\\begin{itemize}
\\item A fortran compiler.
\\item An implementation of BLAS/LAPACK.
\\end{itemize}
Under Ubuntu (14.04 LTS), the installation of the dependencies can be performed using the following command:
\\begin{command}
> sudo apt-get install gfortran libatlas-base-dev
\\end{command}
")
diff --git a/packages/25_damage_non_local.cmake b/packages/25_damage_non_local.cmake
index a2576d36f..9bd94e95a 100644
--- a/packages/25_damage_non_local.cmake
+++ b/packages/25_damage_non_local.cmake
@@ -1,69 +1,63 @@
#===============================================================================
# @file damage_non_local.cmake
#
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Fri Jun 15 13:48:37 2012
#
# @brief package description for non-local materials
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
option(AKANTU_DAMAGE_NON_LOCAL "Package for Non-local damage constitutives laws Akantu" OFF)
-add_internal_package_dependencies(damage_non_local extra_materials)
add_external_package_dependencies(damage_non_local lapack)
set(AKANTU_DAMAGE_NON_LOCAL_FILES
- model/solid_mechanics/materials/weight_function.cc
- model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local.cc
- model/solid_mechanics/materials/material_damage/material_mazars_non_local.cc
-
- model/solid_mechanics/materials/material_non_local_includes.hh
- model/solid_mechanics/materials/weight_function.hh
+ model/solid_mechanics/materials/material_damage/material_damage_non_local.hh
model/solid_mechanics/materials/material_damage/material_marigo_non_local.hh
+ model/solid_mechanics/materials/material_damage/material_marigo_non_local_inline_impl.cc
+ model/solid_mechanics/materials/material_damage/material_mazars_non_local.cc
model/solid_mechanics/materials/material_damage/material_mazars_non_local.hh
- model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local.hh
- model/solid_mechanics/materials/material_damage/material_brittle_non_local.hh
- model/solid_mechanics/materials/material_damage/material_damage_iterative_non_local.hh
model/solid_mechanics/materials/material_non_local.hh
+ model/solid_mechanics/materials/material_non_local_includes.hh
model/solid_mechanics/materials/material_non_local_inline_impl.cc
- model/solid_mechanics/materials/weight_function_tmpl.hh
- model/solid_mechanics/materials/material_damage/material_marigo_non_local_inline_impl.cc
- model/solid_mechanics/materials/material_damage/material_damage_non_local.hh
- model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local_inline_impl.cc
- model/solid_mechanics/materials/material_damage/material_brittle_non_local_inline_impl.cc
- model/solid_mechanics/materials/material_damage/material_damage_iterative_non_local_inline_impl.cc
+ model/solid_mechanics/materials/weight_function.cc
+ model/solid_mechanics/materials/weight_function.hh
+ model/solid_mechanics/materials/weight_function_tmpl.hh
synchronizer/grid_synchronizer.cc
synchronizer/grid_synchronizer.hh
)
set(AKANTU_DAMAGE_NON_LOCAL_TESTS
test_material_damage_non_local
test_grid_synchronizer
)
-
+
+set(AKANTU_DAMAGE_NON_LOCAL_MANUAL_FILES
+ manual-appendix-materials-non-local.tex)
+
set(AKANTU_DAMAGE_NON_LOCAL_DOCUMENTATION "
This package activates the non local damage feature of AKANTU
")
diff --git a/packages/25_damage_non_local.cmake b/packages/25_damage_non_local_extra.cmake
similarity index 60%
copy from packages/25_damage_non_local.cmake
copy to packages/25_damage_non_local_extra.cmake
index a2576d36f..48adb2391 100644
--- a/packages/25_damage_non_local.cmake
+++ b/packages/25_damage_non_local_extra.cmake
@@ -1,69 +1,54 @@
#===============================================================================
# @file damage_non_local.cmake
#
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Fri Jun 15 13:48:37 2012
#
# @brief package description for non-local materials
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
-option(AKANTU_DAMAGE_NON_LOCAL "Package for Non-local damage constitutives laws Akantu" OFF)
+option(AKANTU_DAMAGE_NON_LOCAL_EXTRA "Package for Non-local damage constitutives laws Akantu" OFF)
-add_internal_package_dependencies(damage_non_local extra_materials)
-add_external_package_dependencies(damage_non_local lapack)
+add_internal_package_dependencies(damage_non_local_extra extra_materials)
+add_external_package_dependencies(damage_non_local_extra damage_non_local)
-set(AKANTU_DAMAGE_NON_LOCAL_FILES
- model/solid_mechanics/materials/weight_function.cc
+set(AKANTU_DAMAGE_NON_LOCAL_EXTRA_FILES
model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local.cc
- model/solid_mechanics/materials/material_damage/material_mazars_non_local.cc
- model/solid_mechanics/materials/material_non_local_includes.hh
- model/solid_mechanics/materials/weight_function.hh
- model/solid_mechanics/materials/material_damage/material_marigo_non_local.hh
- model/solid_mechanics/materials/material_damage/material_mazars_non_local.hh
model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local.hh
model/solid_mechanics/materials/material_damage/material_brittle_non_local.hh
model/solid_mechanics/materials/material_damage/material_damage_iterative_non_local.hh
- model/solid_mechanics/materials/material_non_local.hh
- model/solid_mechanics/materials/material_non_local_inline_impl.cc
- model/solid_mechanics/materials/weight_function_tmpl.hh
-
- model/solid_mechanics/materials/material_damage/material_marigo_non_local_inline_impl.cc
- model/solid_mechanics/materials/material_damage/material_damage_non_local.hh
model/solid_mechanics/materials/material_damage/material_vreepeerlings_non_local_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_brittle_non_local_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_damage_iterative_non_local_inline_impl.cc
- synchronizer/grid_synchronizer.cc
- synchronizer/grid_synchronizer.hh
+ model/solid_mechanics/materials/material_non_local_extra_includes.hh
)
-set(AKANTU_DAMAGE_NON_LOCAL_TESTS
- test_material_damage_non_local
- test_grid_synchronizer
+set(AKANTU_DAMAGE_NON_LOCAL_EXTRA_TESTS
)
-
-set(AKANTU_DAMAGE_NON_LOCAL_DOCUMENTATION "
+
+set(AKANTU_DAMAGE_NON_LOCAL_EXTRA_DOCUMENTATION "
This package activates the non local damage feature of AKANTU
")
diff --git a/packages/30_extra_materials.cmake b/packages/30_extra_materials.cmake
index 6d9ccac36..4cd6a4291 100644
--- a/packages/30_extra_materials.cmake
+++ b/packages/30_extra_materials.cmake
@@ -1,86 +1,87 @@
#===============================================================================
# @file extra_materials.cmake
#
# @author Nicolas Richart <nicolas.richart@epfl.ch>
#
# @date Wed Oct 31 16:24:42 2012
#
# @brief package description for extra materials list
#
# @section LICENSE
#
# Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# Akantu is free software: you can redistribute it and/or modify it under the
# terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
# A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with Akantu. If not, see <http://www.gnu.org/licenses/>.
#
#===============================================================================
option(AKANTU_EXTRA_MATERIALS "Add the extra list of materials in Akantu" OFF)
add_external_package_dependencies(extra_materials lapack)
set(AKANTU_EXTRA_MATERIALS_FILES
model/solid_mechanics/materials/material_extra_includes.hh
model/solid_mechanics/materials/material_damage/material_brittle.cc
model/solid_mechanics/materials/material_damage/material_brittle.hh
model/solid_mechanics/materials/material_damage/material_brittle_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_damage_iterative.cc
model/solid_mechanics/materials/material_damage/material_damage_iterative.hh
model/solid_mechanics/materials/material_damage/material_damage_iterative_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_damage_linear.cc
model/solid_mechanics/materials/material_damage/material_damage_linear.hh
model/solid_mechanics/materials/material_damage/material_damage_linear_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_vreepeerlings.hh
model/solid_mechanics/materials/material_damage/material_vreepeerlings_inline_impl.cc
model/solid_mechanics/materials/material_damage/material_vreepeerlings_tmpl.hh
model/solid_mechanics/materials/material_elastic_linear_anisotropic.cc
model/solid_mechanics/materials/material_elastic_linear_anisotropic.hh
model/solid_mechanics/materials/material_elastic_orthotropic.cc
model/solid_mechanics/materials/material_elastic_orthotropic.hh
model/solid_mechanics/materials/material_plastic/material_viscoplastic.cc
model/solid_mechanics/materials/material_plastic/material_viscoplastic.hh
model/solid_mechanics/materials/material_plastic/material_viscoplastic_inline_impl.cc
model/solid_mechanics/materials/material_viscoelastic/material_stiffness_proportional.cc
model/solid_mechanics/materials/material_viscoelastic/material_stiffness_proportional.hh
)
set(AKANTU_EXTRA_MATERIALS_TESTS
)
set(AKANTU_EXTRA_MATERIALS_DOC
manual/manual-extra_materials.tex
+ manual-appendix-materials-extra-materials.tex
)
set(AKANTU_EXTRA_MATERIALS_MANUAL_FILES
manual-extra_materials.tex
figures/stress_strain_visco.pdf
)
set(AKANTU_EXTRA_MATERIALS_DOCUMENTATION "
This package activates additional constitutive laws:
\\begin{itemize}
\\item Linear anisotropy
\\item Linear orthotropy
\\item Visco-plastic
\\end{itemize}
" )
diff --git a/src/common/aka_blas_lapack.hh b/src/common/aka_blas_lapack.hh
index f5a68884a..8c38738e8 100644
--- a/src/common/aka_blas_lapack.hh
+++ b/src/common/aka_blas_lapack.hh
@@ -1,317 +1,316 @@
/**
* @file aka_blas_lapack.hh
*
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date Wed Aug 04 10:58:42 2010
*
* @brief Interface of the Fortran BLAS/LAPACK libraries
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef __AKANTU_AKA_BLAS_LAPACK_HH__
#define __AKANTU_AKA_BLAS_LAPACK_HH__
/* -------------------------------------------------------------------------- */
#if defined(AKANTU_USE_BLAS) || defined(AKANTU_USE_LAPACK)
# include "aka_fortran_mangling.hh"
#endif //AKANTU_USE_BLAS
#ifdef AKANTU_USE_BLAS
extern "C" {
/* ------------------------------------------------------------------------ */
/* Double precision */
/* ------------------------------------------------------------------------ */
//LEVEL 1
double AKA_FC_GLOBAL(ddot, DDOT)(int *, double *, int *, double *, int *);
+
//LEVEL 2
int AKA_FC_GLOBAL(dgemv, DGEMV)(char *, int *, int *, double *, double *, int *,
double *, int *, double *, double *, int *);
//LEVEL 3
int AKA_FC_GLOBAL(dgemm, DGEMM)(char *, char *, int *, int *, int *, double *,
double *, int *, double *, int *, double *,
double *, int *);
-
/* ------------------------------------------------------------------------ */
/* Simple precision */
/* ------------------------------------------------------------------------ */
//LEVEL 1
float AKA_FC_GLOBAL(sdot, SDOT)(int *, float *, int *, float *, int *);
//LEVEL 2
int AKA_FC_GLOBAL(sgemv, SGEMV)(char *, int *, int *, float *, float *, int *,
float *, int *, float *, float *, int *);
//LEVEL 3
int AKA_FC_GLOBAL(sgemm, SGEMM)(char *, char *, int *, int *, int *, float *,
float *, int *, float *, int *, float *,
float *, int *);
-
}
#endif
__BEGIN_AKANTU__
#if defined(__INTEL_COMPILER)
//#pragma warning ( disable : 383 )
#elif defined (__clang__) // test clang to be sure that when we test for gnu it is only gnu
#elif (defined(__GNUC__) || defined(__GNUG__))
# define GCC_VERSION (__GNUC__ * 10000 + __GNUC_MINOR__ * 100 + __GNUC_PATCHLEVEL__)
# if GCC_VERSION > 40600
# pragma GCC diagnostic push
# endif
# pragma GCC diagnostic ignored "-Wunused"
#endif
template<typename T>
inline T aka_dot(int *n, T *x, int *incx, T *y, int *incy) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated BLAS in the compilation options!");
}
template<typename T>
inline int aka_gemv(char *trans, int *m, int *n, T *
alpha, T *a, int *lda, T *x, int *incx,
T *beta, T *y, int *incy) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated BLAS in the compilation options!");
}
template<typename T>
inline int aka_gemm(char *transa, char *transb,
int *m, int *n, int *k,
T *alpha, T *a, int *lda,
T *b, int *ldb,
T *beta, T *c, int *ldc) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated BLAS in the compilation options!");
}
#if defined(AKANTU_USE_BLAS)
template<>
inline double aka_dot<double>(int *n, double *x, int *incx, double *y, int *incy) {
return AKA_FC_GLOBAL(ddot, DDOT)(n, x, incx, y, incy);
}
template<>
inline int aka_gemv<double>(char *trans, int *m, int *n, double *
alpha, double *a, int *lda, double *x, int *incx,
double *beta, double *y, int *incy) {
return AKA_FC_GLOBAL(dgemv, DGEMV)(trans, m, n, alpha, a, lda, x, incx,
beta, y, incy);
}
template<>
inline int aka_gemm<double>(char *transa, char *transb,
int *m, int *n, int *k,
double *alpha, double *a, int *lda,
double *b, int *ldb,
double *beta, double *c, int *ldc) {
return AKA_FC_GLOBAL(dgemm, DGEMM)(transa, transb, m, n, k, alpha, a, lda,
b, ldb, beta, c, ldc);
}
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
template<>
inline float aka_dot<float>(int *n, float *x, int *incx, float *y, int *incy) {
return AKA_FC_GLOBAL(sdot, SDOT)(n, x, incx, y, incy);
}
template<>
inline int aka_gemv<float>(char *trans, int *m, int *n, float *
alpha, float *a, int *lda, float *x, int *incx,
float *beta, float *y, int *incy) {
return AKA_FC_GLOBAL(sgemv, SGEMV)(trans, m, n, alpha, a, lda, x, incx,
beta, y, incy);
}
template<>
inline int aka_gemm<float>(char *transa, char *transb,
int *m, int *n, int *k,
float *alpha, float *a, int *lda,
float *b, int *ldb,
float *beta, float *c, int *ldc) {
return AKA_FC_GLOBAL(sgemm, SGEMM)(transa, transb, m, n, k, alpha, a, lda,
b, ldb, beta, c, ldc);
}
#endif
__END_AKANTU__
#ifdef AKANTU_USE_LAPACK
extern "C" {
/* ------------------------------------------------------------------------ */
/* Double general matrix */
/* ------------------------------------------------------------------------ */
/// compute the eigenvalues/vectors
void AKA_FC_GLOBAL(dgeev, DGEEV)(char* jobvl, char* jobvr, int* n, double* a,
int* lda, double* wr, double* wi, double* vl, int* ldvl,
double* vr, int* ldvr, double* work, int* lwork, int* info);
/// LU decomposition of a general matrix
void AKA_FC_GLOBAL(dgetrf, DGETRF)(int* m, int *n,
double* a, int* lda,
int* ipiv, int* info);
/// generate inverse of a matrix given its LU decomposition
void AKA_FC_GLOBAL(dgetri, DGETRI)(int* n, double* a, int* lda,
int* ipiv, double* work, int* lwork, int* info);
/// solving A x = b using a LU factorization
void AKA_FC_GLOBAL(dgetrs, DGETRS)(char * trans, int * n, int * nrhs,
double * A, int * lda, int * ipiv,
double * b, int * ldb, int * info);
/* ------------------------------------------------------------------------ */
/* Simple general matrix */
/* ------------------------------------------------------------------------ */
/// compute the eigenvalues/vectors
void AKA_FC_GLOBAL(sgeev, SGEEV)(char* jobvl, char* jobvr, int* n, float* a,
int* lda, float* wr, float* wi, float* vl, int* ldvl,
float* vr, int* ldvr, float* work, int* lwork, int* info);
/// LU decomposition of a general matrix
void AKA_FC_GLOBAL(sgetrf, SGETRF)(int* m, int *n,
float* a, int* lda,
int* ipiv, int* info);
/// generate inverse of a matrix given its LU decomposition
void AKA_FC_GLOBAL(sgetri, SGETRI)(int* n, float* a, int* lda,
int* ipiv, float* work, int* lwork, int* info);
/// solving A x = b using a LU factorization
void AKA_FC_GLOBAL(sgetrs, SGETRS)(char * trans, int * n, int * nrhs,
float * A, int * lda, int * ipiv,
float * b, int * ldb, int * info);
}
#endif //AKANTU_USE_LAPACK
__BEGIN_AKANTU__
template<typename T>
inline void aka_geev(char* jobvl, char* jobvr, int* n, T* a,
int* lda, T* wr, T* wi, T* vl, int* ldvl,
T* vr, int* ldvr, T* work, int* lwork, int* info) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated LAPACK in the compilation options!");
}
template<typename T>
inline void aka_getrf(int* m, int *n,
T* a, int* lda,
int* ipiv, int* info) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated LAPACK in the compilation options!");
}
template<typename T>
inline void aka_getri(int* n, T* a, int* lda,
int* ipiv, T* work, int* lwork, int* info) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated LAPACK in the compilation options!");
}
template<typename T>
inline void aka_getrs(char *trans, int * n, int * nrhs,
T * A, int * lda, int * ipiv,
T * b, int * ldb, int * info) {
AKANTU_DEBUG_ERROR(debug::demangle(typeid(T).name()) << "is not a type recognized, or you didn't activated LAPACK in the compilation options!");
}
#if defined(__INTEL_COMPILER)
//#pragma warning ( disable : 383 )
#elif defined (__clang__) // test clang to be sure that when we test for gnu it is only gnu
#elif defined(__GNUG__)
# if GCC_VERSION > 40600
# pragma GCC diagnostic pop
# else
# pragma GCC diagnostic warning "-Wunused"
# endif
#endif
#ifdef AKANTU_USE_LAPACK
template<>
inline void aka_geev<double>(char* jobvl, char* jobvr, int* n, double* a,
int* lda, double* wr, double* wi, double* vl, int* ldvl,
double* vr, int* ldvr, double* work, int* lwork, int* info) {
AKA_FC_GLOBAL(dgeev, DGEEV)(jobvl, jobvr, n, a,
lda, wr, wi, vl, ldvl,
vr, ldvr, work, lwork, info);
}
template<>
inline void aka_getrf<double>(int* m, int *n,
double* a, int* lda,
int* ipiv, int* info) {
AKA_FC_GLOBAL(dgetrf, DGETRF)(m, n, a, lda, ipiv, info);
}
template<>
inline void aka_getri<double>(int* n, double* a, int* lda,
int* ipiv, double* work, int* lwork, int* info) {
AKA_FC_GLOBAL(dgetri, DGETRI)(n, a, lda, ipiv, work, lwork, info);
}
template<>
inline void aka_getrs<double>(char *trans, int * n, int * nrhs,
double * A, int * lda, int * ipiv,
double * b, int * ldb, int * info) {
AKA_FC_GLOBAL(dgetrs, DGETRS)(trans, n, nrhs, A, lda, ipiv, b, ldb, info);
}
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
template<>
inline void aka_geev<float>(char* jobvl, char* jobvr, int* n, float* a,
int* lda, float* wr, float* wi, float* vl, int* ldvl,
float* vr, int* ldvr, float* work, int* lwork, int* info) {
AKA_FC_GLOBAL(sgeev, SGEEV)(jobvl, jobvr, n, a,
lda, wr, wi, vl, ldvl,
vr, ldvr, work, lwork, info);
}
template<>
inline void aka_getrf<float>(int* m, int *n,
float* a, int* lda,
int* ipiv, int* info) {
AKA_FC_GLOBAL(sgetrf, SGETRF)(m, n, a, lda, ipiv, info);
}
template<>
inline void aka_getri<float>(int* n, float* a, int* lda,
int* ipiv, float* work, int* lwork, int* info) {
AKA_FC_GLOBAL(sgetri, SGETRI)(n, a, lda, ipiv, work, lwork, info);
}
template<>
inline void aka_getrs<float>(char *trans, int * n, int * nrhs,
float * A, int * lda, int * ipiv,
float * b, int * ldb, int * info) {
AKA_FC_GLOBAL(sgetrs, SGETRS)(trans, n, nrhs, A, lda, ipiv, b, ldb, info);
}
#endif
__END_AKANTU__
#endif /* __AKANTU_AKA_BLAS_LAPACK_HH__ */
diff --git a/src/fem/cohesive_element.cc b/src/fe_engine/cohesive_element.cc
similarity index 100%
rename from src/fem/cohesive_element.cc
rename to src/fe_engine/cohesive_element.cc
diff --git a/src/fem/cohesive_element.hh b/src/fe_engine/cohesive_element.hh
similarity index 100%
rename from src/fem/cohesive_element.hh
rename to src/fe_engine/cohesive_element.hh
diff --git a/src/fem/element_class.cc b/src/fe_engine/element_class.cc
similarity index 100%
rename from src/fem/element_class.cc
rename to src/fe_engine/element_class.cc
diff --git a/src/fem/element_class.hh b/src/fe_engine/element_class.hh
similarity index 100%
rename from src/fem/element_class.hh
rename to src/fe_engine/element_class.hh
diff --git a/src/fem/element_class_igfem.cc b/src/fe_engine/element_class_igfem.cc
similarity index 100%
rename from src/fem/element_class_igfem.cc
rename to src/fe_engine/element_class_igfem.cc
diff --git a/src/fem/element_class_igfem.hh b/src/fe_engine/element_class_igfem.hh
similarity index 100%
rename from src/fem/element_class_igfem.hh
rename to src/fe_engine/element_class_igfem.hh
diff --git a/src/fem/element_class_structural.hh b/src/fe_engine/element_class_structural.hh
similarity index 100%
rename from src/fem/element_class_structural.hh
rename to src/fe_engine/element_class_structural.hh
diff --git a/src/fem/element_class_tmpl.hh b/src/fe_engine/element_class_tmpl.hh
similarity index 100%
rename from src/fem/element_class_tmpl.hh
rename to src/fe_engine/element_class_tmpl.hh
diff --git a/src/fem/element_classes/element_class_bernoulli_beam_inline_impl.cc b/src/fe_engine/element_classes/element_class_bernoulli_beam_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_bernoulli_beam_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_bernoulli_beam_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_hexahedron_8_inline_impl.cc b/src/fe_engine/element_classes/element_class_hexahedron_8_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_hexahedron_8_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_hexahedron_8_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_igfem_triangle_3_inline_impl.cc b/src/fe_engine/element_classes/element_class_igfem_triangle_3_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_igfem_triangle_3_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_igfem_triangle_3_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_pentahedron_6_inline_impl.cc b/src/fe_engine/element_classes/element_class_pentahedron_6_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_pentahedron_6_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_pentahedron_6_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_point_1_inline_impl.cc b/src/fe_engine/element_classes/element_class_point_1_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_point_1_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_point_1_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_quadrangle_4_inline_impl.cc b/src/fe_engine/element_classes/element_class_quadrangle_4_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_quadrangle_4_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_quadrangle_4_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_quadrangle_8_inline_impl.cc b/src/fe_engine/element_classes/element_class_quadrangle_8_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_quadrangle_8_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_quadrangle_8_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_segment_2_inline_impl.cc b/src/fe_engine/element_classes/element_class_segment_2_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_segment_2_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_segment_2_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_segment_3_inline_impl.cc b/src/fe_engine/element_classes/element_class_segment_3_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_segment_3_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_segment_3_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_tetrahedron_10_inline_impl.cc b/src/fe_engine/element_classes/element_class_tetrahedron_10_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_tetrahedron_10_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_tetrahedron_10_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_tetrahedron_4_inline_impl.cc b/src/fe_engine/element_classes/element_class_tetrahedron_4_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_tetrahedron_4_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_tetrahedron_4_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_triangle_3_inline_impl.cc b/src/fe_engine/element_classes/element_class_triangle_3_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_triangle_3_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_triangle_3_inline_impl.cc
diff --git a/src/fem/element_classes/element_class_triangle_6_inline_impl.cc b/src/fe_engine/element_classes/element_class_triangle_6_inline_impl.cc
similarity index 100%
rename from src/fem/element_classes/element_class_triangle_6_inline_impl.cc
rename to src/fe_engine/element_classes/element_class_triangle_6_inline_impl.cc
diff --git a/src/fem/fe_engine.cc b/src/fe_engine/fe_engine.cc
similarity index 100%
rename from src/fem/fe_engine.cc
rename to src/fe_engine/fe_engine.cc
diff --git a/src/fem/fe_engine.hh b/src/fe_engine/fe_engine.hh
similarity index 100%
rename from src/fem/fe_engine.hh
rename to src/fe_engine/fe_engine.hh
diff --git a/src/fem/fe_engine_inline_impl.cc b/src/fe_engine/fe_engine_inline_impl.cc
similarity index 100%
rename from src/fem/fe_engine_inline_impl.cc
rename to src/fe_engine/fe_engine_inline_impl.cc
diff --git a/src/fem/fe_engine_template.hh b/src/fe_engine/fe_engine_template.hh
similarity index 100%
rename from src/fem/fe_engine_template.hh
rename to src/fe_engine/fe_engine_template.hh
diff --git a/src/fem/fe_engine_template_cohesive.cc b/src/fe_engine/fe_engine_template_cohesive.cc
similarity index 100%
rename from src/fem/fe_engine_template_cohesive.cc
rename to src/fe_engine/fe_engine_template_cohesive.cc
diff --git a/src/fem/fe_engine_template_tmpl.hh b/src/fe_engine/fe_engine_template_tmpl.hh
similarity index 100%
rename from src/fem/fe_engine_template_tmpl.hh
rename to src/fe_engine/fe_engine_template_tmpl.hh
diff --git a/src/fem/geometrical_data_tmpl.hh b/src/fe_engine/geometrical_data_tmpl.hh
similarity index 100%
rename from src/fem/geometrical_data_tmpl.hh
rename to src/fe_engine/geometrical_data_tmpl.hh
diff --git a/src/fem/geometrical_element.cc b/src/fe_engine/geometrical_element.cc
similarity index 100%
rename from src/fem/geometrical_element.cc
rename to src/fe_engine/geometrical_element.cc
diff --git a/src/fem/igfem_element.cc b/src/fe_engine/igfem_element.cc
similarity index 100%
rename from src/fem/igfem_element.cc
rename to src/fe_engine/igfem_element.cc
diff --git a/src/fem/igfem_element.hh b/src/fe_engine/igfem_element.hh
similarity index 100%
rename from src/fem/igfem_element.hh
rename to src/fe_engine/igfem_element.hh
diff --git a/src/fem/integration_element.cc b/src/fe_engine/integration_element.cc
similarity index 100%
rename from src/fem/integration_element.cc
rename to src/fe_engine/integration_element.cc
diff --git a/src/fem/integrator.hh b/src/fe_engine/integrator.hh
similarity index 100%
rename from src/fem/integrator.hh
rename to src/fe_engine/integrator.hh
diff --git a/src/fem/integrator_gauss.hh b/src/fe_engine/integrator_gauss.hh
similarity index 100%
rename from src/fem/integrator_gauss.hh
rename to src/fe_engine/integrator_gauss.hh
diff --git a/src/fem/integrator_gauss_inline_impl.cc b/src/fe_engine/integrator_gauss_inline_impl.cc
similarity index 100%
rename from src/fem/integrator_gauss_inline_impl.cc
rename to src/fe_engine/integrator_gauss_inline_impl.cc
diff --git a/src/fem/interpolation_element.cc b/src/fe_engine/interpolation_element.cc
similarity index 100%
rename from src/fem/interpolation_element.cc
rename to src/fe_engine/interpolation_element.cc
diff --git a/src/fem/interpolation_element_tmpl.hh b/src/fe_engine/interpolation_element_tmpl.hh
similarity index 100%
rename from src/fem/interpolation_element_tmpl.hh
rename to src/fe_engine/interpolation_element_tmpl.hh
diff --git a/src/fem/shape_cohesive.hh b/src/fe_engine/shape_cohesive.hh
similarity index 100%
rename from src/fem/shape_cohesive.hh
rename to src/fe_engine/shape_cohesive.hh
diff --git a/src/fem/shape_cohesive_inline_impl.cc b/src/fe_engine/shape_cohesive_inline_impl.cc
similarity index 100%
rename from src/fem/shape_cohesive_inline_impl.cc
rename to src/fe_engine/shape_cohesive_inline_impl.cc
diff --git a/src/fem/shape_functions.hh b/src/fe_engine/shape_functions.hh
similarity index 100%
rename from src/fem/shape_functions.hh
rename to src/fe_engine/shape_functions.hh
diff --git a/src/fem/shape_functions_inline_impl.cc b/src/fe_engine/shape_functions_inline_impl.cc
similarity index 100%
rename from src/fem/shape_functions_inline_impl.cc
rename to src/fe_engine/shape_functions_inline_impl.cc
diff --git a/src/fem/shape_igfem.hh b/src/fe_engine/shape_igfem.hh
similarity index 100%
rename from src/fem/shape_igfem.hh
rename to src/fe_engine/shape_igfem.hh
diff --git a/src/fem/shape_igfem_inline_impl.cc b/src/fe_engine/shape_igfem_inline_impl.cc
similarity index 100%
rename from src/fem/shape_igfem_inline_impl.cc
rename to src/fe_engine/shape_igfem_inline_impl.cc
diff --git a/src/fem/shape_lagrange.cc b/src/fe_engine/shape_lagrange.cc
similarity index 100%
rename from src/fem/shape_lagrange.cc
rename to src/fe_engine/shape_lagrange.cc
diff --git a/src/fem/shape_lagrange.hh b/src/fe_engine/shape_lagrange.hh
similarity index 100%
rename from src/fem/shape_lagrange.hh
rename to src/fe_engine/shape_lagrange.hh
diff --git a/src/fem/shape_lagrange_inline_impl.cc b/src/fe_engine/shape_lagrange_inline_impl.cc
similarity index 100%
rename from src/fem/shape_lagrange_inline_impl.cc
rename to src/fe_engine/shape_lagrange_inline_impl.cc
diff --git a/src/fem/shape_linked.cc b/src/fe_engine/shape_linked.cc
similarity index 100%
rename from src/fem/shape_linked.cc
rename to src/fe_engine/shape_linked.cc
diff --git a/src/fem/shape_linked.hh b/src/fe_engine/shape_linked.hh
similarity index 100%
rename from src/fem/shape_linked.hh
rename to src/fe_engine/shape_linked.hh
diff --git a/src/fem/shape_linked_inline_impl.cc b/src/fe_engine/shape_linked_inline_impl.cc
similarity index 100%
rename from src/fem/shape_linked_inline_impl.cc
rename to src/fe_engine/shape_linked_inline_impl.cc
diff --git a/src/model/solid_mechanics/material_list.hh b/src/model/solid_mechanics/material_list.hh
index 31318ce5e..09f34509c 100644
--- a/src/model/solid_mechanics/material_list.hh
+++ b/src/model/solid_mechanics/material_list.hh
@@ -1,73 +1,80 @@
/**
* @file material_list.hh
*
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date Tue Oct 29 10:08:25 2013
*
* @brief List of materials and all includes
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef __AKANTU_MATERIAL_LIST_HH__
#define __AKANTU_MATERIAL_LIST_HH__
#include "aka_config.hh"
/* -------------------------------------------------------------------------- */
/* Material list */
/* -------------------------------------------------------------------------- */
#ifndef AKANTU_CMAKE_LIST_MATERIALS
#include "material_elastic.hh"
#endif
#define AKANTU_CORE_MATERIAL_LIST \
((2, (elastic , MaterialElastic )))
#if defined(AKANTU_EXTRA_MATERIALS)
# include "material_extra_includes.hh"
#else
# define AKANTU_EXTRA_MATERIAL_LIST
#endif
#if defined(AKANTU_COHESIVE_ELEMENT)
# include "material_cohesive_includes.hh"
#else
# define AKANTU_COHESIVE_MATERIAL_LIST
#endif
#if defined(AKANTU_DAMAGE_NON_LOCAL)
# include "material_non_local_includes.hh"
#else
# define AKANTU_DAMAGE_NON_LOCAL_MATERIAL_LIST
#endif
+#if defined(AKANTU_DAMAGE_NON_LOCAL_EXTRA)
+# include "material_non_local_extra_includes.hh"
+#else
+# define AKANTU_DAMAGE_NON_LOCAL_EXTRA_MATERIAL_LIST
+#endif
+
#define AKANTU_MATERIAL_LIST \
AKANTU_CORE_MATERIAL_LIST \
AKANTU_EXTRA_MATERIAL_LIST \
AKANTU_COHESIVE_MATERIAL_LIST \
- AKANTU_DAMAGE_NON_LOCAL_MATERIAL_LIST
+ AKANTU_DAMAGE_NON_LOCAL_MATERIAL_LIST \
+ AKANTU_DAMAGE_NON_LOCAL_EXTRA_MATERIAL_LIST
#endif /* __AKANTU_MATERIAL_LIST_HH__ */
diff --git a/src/model/solid_mechanics/materials/material_non_local_includes.hh b/src/model/solid_mechanics/materials/material_non_local_extra_includes.hh
similarity index 87%
copy from src/model/solid_mechanics/materials/material_non_local_includes.hh
copy to src/model/solid_mechanics/materials/material_non_local_extra_includes.hh
index fa0beea9e..575124934 100644
--- a/src/model/solid_mechanics/materials/material_non_local_includes.hh
+++ b/src/model/solid_mechanics/materials/material_non_local_extra_includes.hh
@@ -1,64 +1,59 @@
/**
- * @file material_non_local_includes.hh
+ * @file material_non_local_extra_includes.hh
*
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date Wed Oct 31 16:24:42 2012
*
* @brief Non local materials includes
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef AKANTU_CMAKE_LIST_MATERIALS
-# include "material_marigo_non_local.hh"
-# include "material_mazars_non_local.hh"
# include "material_vreepeerlings_non_local.hh"
# include "material_brittle_non_local.hh"
# include "material_damage_iterative_non_local.hh"
#endif
#define AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST \
((stress_wf, (StressBasedWeightFunction ))) \
((damage_wf, (DamagedWeightFunction ))) \
((remove_wf, (RemoveDamagedWeightFunction))) \
((base_wf, (BaseWeightFunction )))
#define AKANTU_POSSIBLE_DAMAGE_PARENT_MATERIALS \
((sls_deviatoric , (RemoveDamagedWeightFunction)(MaterialStandardLinearSolidDeviatoric))) \
((neohookean_base_wf , (BaseWeightFunction)(MaterialNeohookean ))) \
((neohookean_remove_wf, (RemoveDamagedWeightFunction)(MaterialNeohookean ))) \
((elastic , (RemoveDamagedWeightFunction)(MaterialElastic )))
#define AKANTU_MATERIAL_VREEPEERLINGS_WEIGHT_FUNCTION_TMPL_LIST \
AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST \
((removed_damrate_wf, RemoveDamagedWithDamageRateWeightFunction))
-#define AKANTU_DAMAGE_NON_LOCAL_MATERIAL_LIST \
- ((3, (marigo_non_local , MaterialMarigoNonLocal, \
- AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST))) \
- ((2, (mazars_non_local , MaterialMazarsNonLocal ))) \
+#define AKANTU_DAMAGE_NON_LOCAL_MATERIAL_EXTRA_LIST \
((3, (vreepeerlings_non_local, MaterialVreePeerlingsNonLocal, \
AKANTU_POSSIBLE_DAMAGE_PARENT_MATERIALS))) \
((3, (brittle_non_local , MaterialBrittleNonLocal, \
AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST))) \
((3, (damage_iterative_non_local , MaterialDamageIterativeNonLocal, \
AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST)))
diff --git a/src/model/solid_mechanics/materials/material_non_local_includes.hh b/src/model/solid_mechanics/materials/material_non_local_includes.hh
index fa0beea9e..a967e6b47 100644
--- a/src/model/solid_mechanics/materials/material_non_local_includes.hh
+++ b/src/model/solid_mechanics/materials/material_non_local_includes.hh
@@ -1,64 +1,45 @@
/**
* @file material_non_local_includes.hh
*
* @author Nicolas Richart <nicolas.richart@epfl.ch>
*
* @date Wed Oct 31 16:24:42 2012
*
* @brief Non local materials includes
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef AKANTU_CMAKE_LIST_MATERIALS
# include "material_marigo_non_local.hh"
# include "material_mazars_non_local.hh"
-# include "material_vreepeerlings_non_local.hh"
-# include "material_brittle_non_local.hh"
-# include "material_damage_iterative_non_local.hh"
#endif
#define AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST \
((stress_wf, (StressBasedWeightFunction ))) \
((damage_wf, (DamagedWeightFunction ))) \
((remove_wf, (RemoveDamagedWeightFunction))) \
((base_wf, (BaseWeightFunction )))
-#define AKANTU_POSSIBLE_DAMAGE_PARENT_MATERIALS \
- ((sls_deviatoric , (RemoveDamagedWeightFunction)(MaterialStandardLinearSolidDeviatoric))) \
- ((neohookean_base_wf , (BaseWeightFunction)(MaterialNeohookean ))) \
- ((neohookean_remove_wf, (RemoveDamagedWeightFunction)(MaterialNeohookean ))) \
- ((elastic , (RemoveDamagedWeightFunction)(MaterialElastic )))
-
-#define AKANTU_MATERIAL_VREEPEERLINGS_WEIGHT_FUNCTION_TMPL_LIST \
- AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST \
- ((removed_damrate_wf, RemoveDamagedWithDamageRateWeightFunction))
-
#define AKANTU_DAMAGE_NON_LOCAL_MATERIAL_LIST \
((3, (marigo_non_local , MaterialMarigoNonLocal, \
AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST))) \
- ((2, (mazars_non_local , MaterialMazarsNonLocal ))) \
- ((3, (vreepeerlings_non_local, MaterialVreePeerlingsNonLocal, \
- AKANTU_POSSIBLE_DAMAGE_PARENT_MATERIALS))) \
- ((3, (brittle_non_local , MaterialBrittleNonLocal, \
- AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST))) \
- ((3, (damage_iterative_non_local , MaterialDamageIterativeNonLocal, \
- AKANTU_MATERIAL_WEIGHT_FUNCTION_TMPL_LIST)))
+ ((2, (mazars_non_local , MaterialMazarsNonLocal )))
diff --git a/test/test_contact/hertz_2D.cc b/test/test_contact/hertz_2D.cc
index 63999659f..0cbcf9a16 100644
--- a/test/test_contact/hertz_2D.cc
+++ b/test/test_contact/hertz_2D.cc
@@ -1,164 +1,164 @@
/**
* @file hertz_2D.cc
*
* @author Alejandro M. Aragón <alejandro.aragon@epfl.ch>
*
* @date Mon Mar 24 11:30:00 2014
*
* @brief This file tests for the Hertz solution in 2D
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
#include "implicit_contact_manager.hh"
using namespace akantu;
using std::cout;
using std::endl;
using std::setw;
using std::setprecision;
int main(int argc, char *argv[]) {
// set dimension
static const UInt dim = 2;
// type definitions
typedef Point<dim> point_type;
typedef BoundingBox<dim> bbox_type;
typedef SolidMechanicsModel model_type;
typedef ModelElement<model_type> master_type;
typedef ContactData<dim,model_type> contact_type;
initialize("steel.dat", argc, argv);
// create mesh
Mesh mesh(dim);
// read mesh
mesh.read("hertz_2D.msh");
// create model
model_type model(mesh);
SolidMechanicsModelOptions opt(_static);
// initialize material
model.initFull(opt);
-
+ model.updateCurrentPosition();
// create data structure that holds contact data
contact_type cd(argc, argv, model);
// optimal value of penalty multiplier
cd[Alpha] = 0.4;
// use bounding box to minimize slave-master pairs
Real r0 = 0.5;
Real r1 = 0.15;
point_type c1(-r0/2,-r1/2);
point_type c2( r0/2, r1/2);
bbox_type bb(c1, c2);
// get physical names from mesh
Array<Real> &coords = mesh.getNodes();
mesh.createGroupsFromMeshData<std::string>("physical_names");
// compute areas for slave nodes that are used for the computation of contact pressures
model.applyBC(BC::Neumann::FromHigherDim(Matrix<Real>::eye(2,1.)), "contact_surface");
// NOTE: the areas are computed by assigning a unit pressure to the contact surface,
// then the magnitude of the resulting force vector at nodes gives its associated area
Array<Real>& areas = model.getForce();
// add slave-master pairs and store slave node areas
ElementGroup &eg = mesh.getElementGroup("contact_surface");
ElementGroup &rs = mesh.getElementGroup("rigid");
for (auto nit = eg.node_begin(); nit != eg.node_end(); ++nit) {
// get point of slave node
point_type n(&coords(*nit));
// process only if within bounding box
if (bb & n) {
// loop over element types
for (ElementGroup::type_iterator tit = rs.firstType(); tit != rs.lastType(); ++tit)
// loop over elements of the rigid surface
for (ElementGroup::const_element_iterator it = rs.element_begin(*tit);
it != rs.element_end(*tit); ++it) {
// create master element
master_type m(model, _segment_2, *it);
assert(has_projection(n,m.point<2>(0),m.point<2>(1)));
// add slave-master pair
cd.addPair(*nit,m);
}
// compute and add area to slave node
Real a = 0.;
for (UInt i=0; i<dim; ++i)
a += pow(areas(*nit, i),2.);
cd.addArea(*nit, sqrt(a));
}
}
// clear force vector before the start of the simulation
areas.clear();
// output contact data info
cout<<cd;
// apply boundary conditions
model.applyBC(BC::Dirichlet::FixedValue(0., BC::_x), "rigid");
model.applyBC(BC::Dirichlet::FixedValue(0., BC::_y), "rigid");
model.getBlockedDOFs()(7,0) = true;
Real data[3][50]; // store results for printing
Real step = 0.001; // top displacement increment
Real Delta = 0.05; // maximum imposed displacement
size_t k=0;
// loop over displacement increments
for (Real delta = step; delta <= Delta+step; delta += step) {
// apply displacement at the top
model.applyBC(BC::Dirichlet::FixedValue(-delta, BC::_y), "top");
// solve contact step (no need to call solve on the model object)
cd.solveContactStep();
data[0][k] = delta;
data[1][k] = cd.getForce();
data[2][k] = cd.getMaxPressure();
++k;
}
// print results
size_t w = 10;
cout<<setprecision(2);
cout<<setw(w)<<"\nDisp."<<setw(w)<<"Force"<<setw(w)<<"Max pressure"<<endl;
for (int i=0; i<50; ++i)
cout<<setw(w)<<data[0][i]<<setw(w)<<data[1][i]<<setw(w)<<data[2][i]<<endl;
// finalize simulation
finalize();
return EXIT_SUCCESS;
}
diff --git a/test/test_contact/hertz_3D.cc b/test/test_contact/hertz_3D.cc
index 9c1d6d745..645332352 100644
--- a/test/test_contact/hertz_3D.cc
+++ b/test/test_contact/hertz_3D.cc
@@ -1,209 +1,210 @@
/**
* @file hertz_3D.cc
*
* @author Alejandro M. Aragón <alejandro.aragon@epfl.ch>
*
* @date Fri Apr 4 11:30:00 2014
*
* @brief This file tests for the Hertz solution in 3D
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* Akantu is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Akantu is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Akantu. If not, see <http://www.gnu.org/licenses/>.
*
*/
#include "implicit_contact_manager.hh"
using namespace akantu;
using std::cout;
using std::endl;
using std::setw;
using std::setprecision;
struct Assignator : public ContactData<3,SolidMechanicsModel>::SearchBase {
- typedef Point <3> point_type;
- typedef SolidMechanicsModel model_type;
- typedef ModelElement <model_type> master_type;
+ typedef Point <3> point_type;
+ typedef SolidMechanicsModel model_type;
+ typedef ModelElement <model_type> master_type;
- model_type &model_;
+ model_type &model_;
- Assignator(model_type& model) : model_(model) {
- }
+ Assignator(model_type& model) : model_(model) {
+ }
- virtual int search(const Real *ptr) {
- point_type p(ptr);
+ virtual int search(const Real *ptr) {
+ point_type p(ptr);
- ElementGroup &rs = model_.getMesh().getElementGroup("rigid_surface");
+ ElementGroup &rs = model_.getMesh().getElementGroup("rigid_surface");
- for (ElementGroup::type_iterator tit = rs.firstType(); tit != rs.lastType(); ++tit)
- for (ElementGroup::const_element_iterator it = rs.element_begin(*tit);
- it != rs.element_end(*tit); ++it) {
- master_type m(model_, _triangle_3, *it);
- if (point_has_projection_to_triangle(p, m.point <3>(0), m.point <3>(1), m.point <3>(2)))
- return m.element;
- }
- return -1;
- }
+ for (ElementGroup::type_iterator tit = rs.firstType(); tit != rs.lastType(); ++tit)
+ for (ElementGroup::const_element_iterator it = rs.element_begin(*tit);
+ it != rs.element_end(*tit); ++it) {
+ master_type m(model_, _triangle_3, *it);
+ if (point_has_projection_to_triangle(p, m.point <3>(0), m.point <3>(1), m.point <3>(2)))
+ return m.element;
+ }
+ return -1;
+ }
};
int main(int argc, char *argv[]) {
- // set dimension
- static const UInt dim = 3;
+ // set dimension
+ static const UInt dim = 3;
- // type definitions
- typedef Point <dim> point_type;
- typedef BoundingBox <dim> bbox_type;
- typedef SolidMechanicsModel model_type;
- typedef ModelElement <model_type> master_type;
- typedef ContactData <dim, model_type> contact_type;
+ // type definitions
+ typedef Point <dim> point_type;
+ typedef BoundingBox <dim> bbox_type;
+ typedef SolidMechanicsModel model_type;
+ typedef ModelElement <model_type> master_type;
+ typedef ContactData <dim, model_type> contact_type;
- initialize("steel.dat", argc, argv);
+ initialize("steel.dat", argc, argv);
- // create mesh
- Mesh mesh(dim);
+ // create mesh
+ Mesh mesh(dim);
- // read mesh
- mesh.read("hertz_3D.msh");
+ // read mesh
+ mesh.read("hertz_3D.msh");
- // create model
- model_type model(mesh);
- SolidMechanicsModelOptions opt(_static);
+ // create model
+ model_type model(mesh);
+ SolidMechanicsModelOptions opt(_static);
- // initialize material
- model.initFull(opt);
+ // initialize material
+ model.initFull(opt);
+ model.updateCurrentPosition();
- // create data structure that holds contact data
- contact_type cd(argc, argv, model);
+ // create data structure that holds contact data
+ contact_type cd(argc, argv, model);
- // optimal value of penalty multiplier
- cd[Alpha] = 0.05;
- cd[Uzawa_tol] = 1.e-2;
- cd[Newton_tol] = 1.e-2;
+ // optimal value of penalty multiplier
+ cd[Alpha] = 0.05;
+ cd[Uzawa_tol] = 1.e-2;
+ cd[Newton_tol] = 1.e-2;
-// // set Paraview output resluts
-// std::string para_file = "contact";
-// model.setBaseName(para_file);
-// model.addDumpFieldVector("displacement");
+ // // set Paraview output resluts
+ // std::string para_file = "contact";
+ // model.setBaseName(para_file);
+ // model.addDumpFieldVector("displacement");
- ///////////////////////////////////////////////////////////////////
-
- // call update current position to be able to call later
- // the function to get current positions
- model.updateCurrentPosition();
-
-
- // get physical names from Gmsh file
- mesh.createGroupsFromMeshData <std::string>("physical_names");
-
- // get areas for the nodes of the circle
- // this is done by applying a unit pressure to the contact surface elements
- model.applyBC(BC::Neumann::FromHigherDim(Matrix <Real>::eye(3, 1.)), "contact_surface");
- Array <Real>& areas = model.getForce();
-
- // loop over contact surface nodes and store node areas
- ElementGroup &eg = mesh.getElementGroup("contact_surface");
- cout << "Adding areas to slave nodes. " << endl;
- for (auto nit = eg.node_begin(); nit != eg.node_end(); ++nit) {
- // compute area contributing to the slave node
- Real a = 0.;
- for (UInt i = 0; i < dim; ++i)
- a += pow(areas(*nit, i), 2.);
- cd.addArea(*nit, sqrt(a));
- }
-
- // set force value to zero
- areas.clear();
-
- // apply boundary conditions for the rigid plane
- model.applyBC(BC::Dirichlet::FixedValue(0., BC::_x), "bottom_body");
- model.applyBC(BC::Dirichlet::FixedValue(0., BC::_y), "bottom_body");
- model.applyBC(BC::Dirichlet::FixedValue(0., BC::_z), "bottom_body");
-
- ElementGroup &cs = mesh.getElementGroup("contact_surface");
-
- Array <Real> &coords = mesh.getNodes();
-
- // set-up bounding box to include slave nodes that lie inside it
- Real l1 = 1.;
- Real l2 = 0.2;
- Real l3 = 1.;
- point_type c1(-l1 / 2, -l2 / 2, -l3 / 2);
- point_type c2(l1 / 2, l2 / 2, l3 / 2);
- bbox_type bb(c1, c2);
-
- // search policy for slave-master pairs
- Assignator a(model);
-
-
- // loop over nodes in contact surface to create contact elements
- cout << "Adding slave-master pairs" << endl;
- for (auto nit = cs.node_begin(); nit != cs.node_end(); ++nit) {
- point_type p(&coords(*nit));
-
- // ignore slave node if it doesn't lie within the bounding box
- if (!(bb & p))
- continue;
-
- int el = a.search(&coords(*nit));
- if (el != -1)
- cd.addPair(*nit, master_type(model, _triangle_3, el));
- }
+ ///////////////////////////////////////////////////////////////////
+
+ // call update current position to be able to call later
+ // the function to get current positions
+ model.updateCurrentPosition();
+
+
+ // get physical names from Gmsh file
+ mesh.createGroupsFromMeshData <std::string>("physical_names");
+
+ // get areas for the nodes of the circle
+ // this is done by applying a unit pressure to the contact surface elements
+ model.applyBC(BC::Neumann::FromHigherDim(Matrix <Real>::eye(3, 1.)), "contact_surface");
+ Array <Real>& areas = model.getForce();
+
+ // loop over contact surface nodes and store node areas
+ ElementGroup &eg = mesh.getElementGroup("contact_surface");
+ cout << "Adding areas to slave nodes. " << endl;
+ for (auto nit = eg.node_begin(); nit != eg.node_end(); ++nit) {
+ // compute area contributing to the slave node
+ Real a = 0.;
+ for (UInt i = 0; i < dim; ++i)
+ a += pow(areas(*nit, i), 2.);
+ cd.addArea(*nit, sqrt(a));
+ }
+
+ // set force value to zero
+ areas.clear();
+
+ // apply boundary conditions for the rigid plane
+ model.applyBC(BC::Dirichlet::FixedValue(0., BC::_x), "bottom_body");
+ model.applyBC(BC::Dirichlet::FixedValue(0., BC::_y), "bottom_body");
+ model.applyBC(BC::Dirichlet::FixedValue(0., BC::_z), "bottom_body");
+
+ ElementGroup &cs = mesh.getElementGroup("contact_surface");
+
+ Array <Real> &coords = mesh.getNodes();
+
+ // set-up bounding box to include slave nodes that lie inside it
+ Real l1 = 1.;
+ Real l2 = 0.2;
+ Real l3 = 1.;
+ point_type c1(-l1 / 2, -l2 / 2, -l3 / 2);
+ point_type c2(l1 / 2, l2 / 2, l3 / 2);
+ bbox_type bb(c1, c2);
+
+ // search policy for slave-master pairs
+ Assignator a(model);
+
+
+ // loop over nodes in contact surface to create contact elements
+ cout << "Adding slave-master pairs" << endl;
+ for (auto nit = cs.node_begin(); nit != cs.node_end(); ++nit) {
+ point_type p(&coords(*nit));
+
+ // ignore slave node if it doesn't lie within the bounding box
+ if (!(bb & p))
+ continue;
+
+ int el = a.search(&coords(*nit));
+ if (el != -1)
+ cd.addPair(*nit, master_type(model, _triangle_3, el));
+ }
-// model.dump();
-
- // block z-disp in extreme points of top surface
- model.getBlockedDOFs()(1, 2) = true;
- model.getBlockedDOFs()(2, 2) = true;
+ // model.dump();
+
+ // block z-disp in extreme points of top surface
+ model.getBlockedDOFs()(1, 2) = true;
+ model.getBlockedDOFs()(2, 2) = true;
- // block x-disp in extreme points of top surface
- model.getBlockedDOFs()(3, 0) = true;
- model.getBlockedDOFs()(4, 0) = true;
-
- const size_t steps = 20;
- Real data[3][steps]; // store results for printing
- Real step = 0.001; // top displacement increment
- size_t k = 0;
-
- for (Real delta = 0; delta <= step * steps; delta += step) {
- // apply displacement to the top surface of the half-sphere
- model.applyBC(BC::Dirichlet::FixedValue(-delta, BC::_y), "top_surface");
-
- // solve contact step, this function also dumps Paraview files
- cd.solveContactStep(&a);
-
- data[0][k] = delta;
- data[1][k] = cd.getForce();
- data[2][k] = cd.getMaxPressure();
- ++k;
- }
-
- // print results
- size_t w = 14;
- cout << setprecision(4);
- cout << endl << setw(w) << "Disp." << setw(w) << "Force" << setw(w) << "Max pressure" << endl;
- for (size_t i = 0; i < steps; ++i)
- cout << setw(w) << data[0][i] << setw(w) << data[1][i] << setw(w) << data[2][i] << endl;
-
- // finalize simulation
- finalize();
- return EXIT_SUCCESS;
+ // block x-disp in extreme points of top surface
+ model.getBlockedDOFs()(3, 0) = true;
+ model.getBlockedDOFs()(4, 0) = true;
+
+ const size_t steps = 20;
+ Real data[3][steps]; // store results for printing
+ Real step = 0.001; // top displacement increment
+ size_t k = 0;
+
+ for (Real delta = 0; delta <= step * steps; delta += step) {
+ // apply displacement to the top surface of the half-sphere
+ model.applyBC(BC::Dirichlet::FixedValue(-delta, BC::_y), "top_surface");
+
+ // solve contact step, this function also dumps Paraview files
+ cd.solveContactStep(&a);
+
+ data[0][k] = delta;
+ data[1][k] = cd.getForce();
+ data[2][k] = cd.getMaxPressure();
+ ++k;
+ }
+
+ // print results
+ size_t w = 14;
+ cout << setprecision(4);
+ cout << endl << setw(w) << "Disp." << setw(w) << "Force" << setw(w) << "Max pressure" << endl;
+ for (size_t i = 0; i < steps; ++i)
+ cout << setw(w) << data[0][i] << setw(w) << data[1][i] << setw(w) << data[2][i] << endl;
+
+ // finalize simulation
+ finalize();
+ return EXIT_SUCCESS;
}