\documentclass[10pt]{article}

\usepackage{fullpage}
\usepackage{setspace}
\usepackage{parskip}
\usepackage{titlesec}
\usepackage[section]{placeins}
\usepackage{xcolor}
\usepackage{breakcites}
\usepackage{lineno}
\usepackage{hyphenat}





\PassOptionsToPackage{hyphens}{url}
\usepackage[colorlinks = true,
            linkcolor = blue,
            urlcolor  = blue,
            citecolor = blue,
            anchorcolor = blue]{hyperref}
\usepackage{etoolbox}
\makeatletter
\patchcmd\@combinedblfloats{\box\@outputbox}{\unvbox\@outputbox}{}{%
  \errmessage{\noexpand\@combinedblfloats could not be patched}%
}%
\makeatother


\usepackage[round]{natbib}
\let\cite\citep




\renewenvironment{abstract}
  {{\bfseries\noindent{\abstractname}\par\nobreak}\footnotesize}
  {\bigskip}

\titlespacing{\section}{0pt}{*3}{*1}
\titlespacing{\subsection}{0pt}{*2}{*0.5}
\titlespacing{\subsubsection}{0pt}{*1.5}{0pt}


\usepackage{authblk}


\usepackage{graphicx}
\usepackage[space]{grffile}
\usepackage{latexsym}
\usepackage{textcomp}
\usepackage{longtable}
\usepackage{tabulary}
\usepackage{booktabs,array,multirow}
\usepackage{amsfonts,amsmath,amssymb}
\providecommand\citet{\cite}
\providecommand\citep{\cite}
\providecommand\citealt{\cite}
% You can conditionalize code for latexml or normal latex using this.
\newif\iflatexml\latexmlfalse
\providecommand{\tightlist}{\setlength{\itemsep}{0pt}\setlength{\parskip}{0pt}}%

\AtBeginDocument{\DeclareGraphicsExtensions{.pdf,.PDF,.eps,.EPS,.png,.PNG,.tif,.TIF,.jpg,.JPG,.jpeg,.JPEG}}

\usepackage[utf8]{inputenc}
\usepackage[ngerman,english]{babel}










\begin{document}

\title{Solving Linear Elasticity Equation on Different Geometry Using FEniCS
Platform}



\author[1]{Hao Lu}%
\affil[1]{University of Memphis}%


\vspace{-1em}



  
  \date{May 27, 2024}


\begingroup
\let\center\flushleft
\let\endcenter\endflushleft
\maketitle
\endgroup








\sloppy


\section*{Abstract}

{\label{785430}}

In this study, linear elasticity equation is solved on different
geometries and under various boundary conditions.~ Classical elasticity
problems like plane stress, uniaxial stress problem are implemented on
regarding geometry. Also, behaviors of relatively complex structures are
computed as well. The key of this project is familarize the work flow
and regarding command using FEniCS to solve partial differencial
equation (PDE). In general, FEniCS is easy to use, and freedom of
controling the computing process is allowed at certain degree through
input script. In this study, some more complex structures suffered from
memory limit of local machine. The ability of FEniCS will be fully
expanded on~high-performance clusters.

\par\null

\section*{Introduction}

{\label{881984}}

FEniCS is an open~source~computing platform for solving partial
differential equations (PDEs). It~enables users to quickly translate
scientific models into efficient finite element
code~\cite{project}.~ Once user specify the the variational form,
domain, boundary condition and the function space, FEniCS autometically
generates the global stiffness matrix and force vector and solve for the
unknown displacement vector. In this study, the linear elasticity
equation is applied to different geometry and boundary condition.~

The main constituents of a finite element method for the solution of a
boundary-value problem are~\cite{BELYTSCHKO_2008} :

\begin{enumerate}
\tightlist
\item
  The variational or weak statement of the problem.
\item
  The approximate solution of the variational equations through the use
  of finite element finctions.
\end{enumerate}

We first start with the methmatical statement of the problem to be
solve. The governing equation for small elastic deformation of a
domain~\(\Omega\)~ in differencial form (or strong form) can be
written as~\cite{solvers} :

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(-\nabla\cdot\sigma\ =\ f\ in\ \Omega\),

\(\sigma\ =\ \lambda tr\left(\varepsilon\right)I\ +2\mu\varepsilon\),

\(\varepsilon\ =\ \frac{1}{2}\left(\nabla u\ +\left(\nabla u\right)^T\right)\),
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

where~\(\sigma\) is the stress tensor,~\(f\) is the
body force per unit volume,~\(\lambda\) and~\(\mu\)
are~Lam\selectlanguage{ngerman}é's elasticity parameters for material in
domain~\(\Omega\),~\(I\) is the identity
temsor,~\(tr\) is the trace operator on a
tensor,~\(\varepsilon\) is the symmetric strain-rate temsor,
and~\(u\) is the displacement~ vector field.~

Combine the second and the third equation, we have:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(\sigma\ =\ \lambda\left(\nabla\cdot u\right)I\ +\ \mu\left(\nabla u\ +\ \left(\nabla u\right)^T\right)\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

Next, we convert the equation into variational form (or weak form):

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(-\int_{_{\Omega}}^{ }\left(\nabla\cdot\sigma\right)\cdot vdx\ =\ \int_{_{\Omega}}^{ }f\cdot vdx\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}

\(-\int_{\Omega}^{ }\left(\nabla\cdot\sigma\right)\cdot vdx\ =\ \) \(\int_{\Omega}^{ }\sigma\ :\ \nabla vdx\ \)~~\(-\int_{\partial\Omega}^{ }\left(\sigma\cdot n\right)\cdot vds\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

The colon operator is the inner product between
tensors,~\(n\) is the outward unit normal at the
boundary.~\(\sigma\cdot n\) represants the value of traction
vector~\(T\) at the boundary and often prescribed as
boundary condition. And the value of the displacement is assumed as
Dirichlet boundary condition.

We have:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(\int_{\Omega}^{ }\sigma\ :\ \nabla vdx\ =\ \int_{\Omega}^{ }f\cdot vdx\ +\ \int_{\partial\Omega_T}^{ }T\cdot vds\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

The variational formulation can be finally summarized as:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(a\left(u,v\right)\ =\ L\left(v\right)\ \ \ \forall v\ \in\ V\),
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

where:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(a\left(u,v\right)\ =\ \int_{\Omega}^{ }\sigma\left(u\right)\ :\ \nabla vdx\),
\end{quote}
\end{quote}

\(\sigma\left(u\right)\ =\ \lambda\left(\nabla\cdot u\right)I\ +\ \mu\left(\nabla u\ +\ \left(\nabla u\right)^T\right)\),

\begin{quote}
\(L\left(v\right)\ =\ \int_{\Omega}^{ }f\cdot vdx\ +\ \int_{\partial\Omega_T}^{ }T\cdot vds\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

Once the variational form is defined, solving for
displacement~\(u\) is relatively striaight forward using
FEniCS. As soon as~\(u\) is computed, we can
compute~\emph{von Mises stress}~defined as:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(\sigma_M=\sqrt{\frac{3}{2}s\ :\ s}\)~
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

where~\(s\) is the deviatoric stress tensor given as:

\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\begin{quote}
\(s\ =\ \sigma-\frac{1}{3}tr\left(\sigma\right)I\)
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}
\end{quote}

The~\emph{von Mises stress}~is often used in determining whether
an~\textbf{isotropic}~meterial~will yield when subjected to a complex
loading condition~\cite{stress}.

In this study, the stress is applied through specify the displacement of
boundaries as Dirichlet boundary conditions. And the problem is scaled
and physical parameters are dimensionsless and on scale of 1.

\subsection*{FEniCS Implementation
Steps}

{\label{765441}}

A typical FEniCS Python code has the following structure:~

\begin{enumerate}
\tightlist
\item
  Import essencial libraries (e.g. fenics, dolfin, mshr)
\item
  Define parameters in the equation
\item
  Define domain, create mesh and define function space
\item
  Define boundary condition
\item
  Define variational problem
\item
  Compute solution~
\item
  Save results for visulization
\end{enumerate}

The actrually Python scripts are not included in the report. Please
follow this~\href{https://osf.io/ak4t3/}{link}~to check the scripts.

\section*{}

{\label{588642}}

\section*{Results}

{\label{588642}}

\subsection*{Uniaxial Stress}

{\label{409000}}

In the state of uniaxial stress, only one of the principal
stresses~\(\sigma_1\) is nonzero,~\(\sigma_2,\ \sigma_3\) are zero. The
geometry used in this case is a 60*30*10 plane with no body force
subjected and a constant displacement specified by Dirichlet boundary
conditions on left and right boundary. The displacement is set to be -5
on left boundary and 5 on the right boundary. The displacement magnitude
decrease from boundary towards the center (Fig. 1).~ The~\emph{von Mises
stress~}(Fig.2) shows yielding would first occur at the conors, however,
the result can be artifical due to the boundary conditions. Apart from
that, the whole domain shows basically equal~~\emph{von Mises
stress,~}which is more convicing under isotropic material
setting.\emph{~}

\par\null\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-05-09-33/Screen-Shot-2017-12-05-at-05-09-33}
\caption{{Magnitude of displacement~
{\label{489581}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-05-09-42/Screen-Shot-2017-12-05-at-05-09-42}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{936195}}%
}}
\end{center}
\end{figure}

\subsection*{Plane Stress}

{\label{398813}}

The state of plane stress exists when there is only one zero component
of principal stress where~\(\sigma_3=0\),~\(\sigma_1\ne0,\ \sigma_2\ne0\).
The~geometry used in this case is a 10*10*2 plane with no body force
subjected. Dirichlet boundary conditions are spcified as fixed
displacements towards the center of the plane to represant plane stress.
The magnitude of the displacement of front, back, left and right
boundary is 5, move toward the center of the domain. Result shows that
strain decreases towards the center of the plane. And the conors are
where the yield will first occor assuming isotropic material, which,
again, might be artifically caused by boundary condition. The rest of
the domain shows basically equal~Magnitude of~\emph{von Mises stress.~}

\par\null\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-02-59-55/Screen-Shot-2017-12-05-at-02-59-55}
\caption{{Magnitude of displacement~
{\label{284436}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-03-00-08/Screen-Shot-2017-12-05-at-03-00-08}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{353954}}%
}}
\end{center}
\end{figure}

\par\null

\subsection*{Strain Computation in Relatively Complex
Structure}

{\label{141727}}

Next, we move to the more complex structures. The first one (Fig.5
Fig.6) is subjected to displacement along x-direction specified
by~Dirichlet boundary conditions (5 on left boundary, -5 on right
boundary). The displacement of predefined boundaries are largest and
decreases towards the middle of the boundary along the other sides.~\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-11-01/Screen-Shot-2017-12-05-at-04-11-01}
\caption{{Magnitude of displacement~
{\label{848842}}%
}}
\end{center}
\end{figure}

The~Magnitude of~\emph{von Mises stress~}is higher on the free
boundaries compare to the~Dirichlet boundaries.~\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-11-11/Screen-Shot-2017-12-05-at-04-11-11}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{733092}}%
}}
\end{center}
\end{figure}

Figure 7, 8, 9, 10 shows the same domain with overall size 20*20*5,
subjected to different boundary conditions. Figure 7 and 8 shows the
domain with displacement applied to the left and right boundaries
(x-direction, 5 on the left and -5 on the right boundary); the domain in
figure 9 and 10 are subjected to displacement along the y-direction (5
on the front and -5 on the back boundary). The magnitude of displacement
was set to be the same in each model. Due to the structure of the
domain, change in direction of the displacement changes where the
yielding will first occur.~\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-13-54/Screen-Shot-2017-12-05-at-04-13-54}
\caption{{Magnitude of displacement~
{\label{312281}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-14-08/Screen-Shot-2017-12-05-at-04-14-08}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{369565}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-19-30/Screen-Shot-2017-12-05-at-04-19-30}
\caption{{Magnitude of displacement~
{\label{705548}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-19-53/Screen-Shot-2017-12-05-at-04-19-53}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{455348}}%
}}
\end{center}
\end{figure}

Model in figure 11 and 12 has the~Dirichlet boundaries on bottom and
top. The bottom boundary is set to be fixed and the top boundary has a
downward displacement with magnitude of 5. In figure 11, the upper part
of the structure has the maximum displacement, yet the strain is
accumulated mostly along the middle where the structure is much
thinner.~ With wider geometry, the upper part is harder to deform.
Figure 12 further shows that the yielding will first occur at the middle
part of the whole structure.\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-23-33/Screen-Shot-2017-12-06-at-00.09.16}
\caption{{Magnitude of displacement~
{\label{616368}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-23-42/Screen-Shot-2017-12-06-at-00.09.35}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{164447}}%
}}
\end{center}
\end{figure}

Figure 13\textasciitilde{}16 shows the comparison~between two structures
under the same boundary condition. The first model of the two has
exactly the same pillars and the second model has different pillars in
terms of width (ranging from 1 to 4). In both models, the bottom
boundary is set to be fixed where no displacement is allowed and top
boundaries have the same amount of displacement. What's interesting is
that with the same amount of vertical strain, the thinner pillar in the
second model did not show a higher tendency of yielding (or at least not
visible under such loading condition).~~\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-50-59/Screen-Shot-2017-12-05-at-04-50-59}
\caption{{Magnitude of displacement~
{\label{558975}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-51-10/Screen-Shot-2017-12-05-at-04-51-10}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{551882}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-05-01-25/Screen-Shot-2017-12-05-at-05-01-25}
\caption{{Magnitude of displacement~
{\label{145687}}%
}}
\end{center}
\end{figure}\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-04-59-56/Screen-Shot-2017-12-05-at-11.51.51}
\caption{{Magnitude of~\emph{von Mises stress}
{\label{991168}}%
}}
\end{center}
\end{figure}

\section*{}

{\label{246231}}

\section*{}

{\label{246231}}

\section*{Discussion~\& Issue}

{\label{246231}}

Running FEniCS on the local machine through Docker seems to suffer from
the~relatively complex structures. There are several more interesting
geometries implemented in this study, but FEniCS ran into ``out of
memory'' issue by only applying the condition that geometry is
subjecting to body force. The same error occurred when applying traction
or Dirichlet boundary condition on their own with the same geometry.
Looks like applying~simple~mathematical~expression~on a~complex~domain
is more expensive to~solve~compares~to solving a~complex~equation~on
a~simple~domain. Also, by default system~limitation of RAM usage on
Docker, even the computing power of local machine cannot be fully
utilized.~

~\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-00-59-48/Screen-Shot-2017-12-05-at-00-59-48}
\caption{{10*10*10 cube with a hollow sphere structure inside
{\label{144336}}%
}}
\end{center}
\end{figure}

~

\par\null\selectlanguage{english}
\begin{figure}[h!]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/Screen-Shot-2017-12-05-at-03-45-33/Screen-Shot-2017-12-05-at-03-45-33}
\caption{{10*10*2 plane with a hollow cone structure at the center
{\label{648625}}%
}}
\end{center}
\end{figure}

\section*{}

{\label{296461}}

\section*{Conclusion}

{\label{296461}}

FEniCS is a user-friendly platform for solving PDE using finite element
method. People with basic FEM knowledge can relatively easy to
understand and start using it to solve PDE and visualize the result.
Especially after hands-on coding the whole process, one can fully
appreciate~the idea of FEniCS in terms of convenience. However, in this
study, it suffers from complex geometry and failed to converge because
of~``out of memory'' issue when running on local machine. What's more,
this study only applied linear elasticity equation without
discretization and integration~over the time. It's behavior under
applying a~much more complex rheological~law and time discretization is
yet to be tested.~

\selectlanguage{english}
\FloatBarrier
\bibliographystyle{plainnat}
\bibliography{bibliography/converted_to_latex.bib%
}

\end{document}

