Homec4science
Diffusion tamaas (features/netcdf-mpi) Browse

tamaas/doc/sphinx/source/performance.rst
343934da4495features/netcdf-mpi

Last Change
Raw File

performance.rst
View Options

Performance
===========
Parallelism
-----------
Tamaas implements shared-memory parallelism using `thrust
<https://github.com/thrust/thrust>`_. The Thrust backend can be controlled with
the following values of the ``backend`` build option:
``omp``
Thrust uses its OpenMP backend. The number of threads is controlled by
OpenMP.
``cpp``
Thurst does not run in threads (i.e. sequential). This is the recommanded
option if running multiple MPI tasks.
``tbb``
Thrust uses its `TBB
<https://en.wikipedia.org/wiki/Threading_Building_Blocks>`_ backend. Note
that this option is not fully supported by Tamaas.
.. tip:: When using the OpenMP or TBB backend, the number of threads can be
manually controlled by the :py:func:`initialize <tamaas._tamaas.initialize>`
function. When OpenMP is select for the backend, the environment variable
``OMP_NUM_THREADS`` can also be used to set the number of threads.
FFTW has its own system for thread-level parallelism, which can be controlled
via the ``fftw_threads`` option:
``none``
FFTW does not use threads.
``threads``
FFTW uses POSIX/Win32 threads for parallelism.
``omp``
FFTW uses OpenMP.
.. note:: As with the Thrust backend, the number of threads for FFTW can be
controlled with :py:func:`initialize <tamaas._tamaas.initialize>`.
Finally, the boolean variable ``use_mpi`` controls wheter Tamaas is compiled
with MPI-parallelism. If yes, Tamaas will be linked against ``libfftw3_mpi``
regardless of the thread model.
Multi-process parallelism
^^^^^^^^^^^^^^^^^^^^^^^^^
Distributed memory parallelism in Tamaas is implemented with `MPI
<https://www.mpi-forum.org/>`_. Due to the bottleneck role of the fast-Fourier
transform in Tamaas' core routines, the data layout of Tamaas is that of `FFTW
<http://fftw.org/fftw3_doc/MPI-Data-Distribution.html#MPI-Data-Distribution>`_.
Tamaas is somewhat affected by limitations of FFTW, and MPI only works on
systems with a 2D boundary, i.e. ``basic_2d``, ``surface_2d`` and ``volume_2d``
model types (which are the most important anyways, since rough contact mechanics
can yield different scaling laws in 1D).
MPI support in Tamaas is still experimental, but the following parts are tested:
- Rough surface generation
- Surface statistics computation
- Elastic normal contact
- Elastic-plastic contact (with :py:class:`DFSANECXXSolver
<tamaas.nonlinear_solvers.DFSANECXXSolver>`).
Adapting existing scripts to work in an MPI context can require some work,
especially if said scripts rely on numpy and scipy for pre- and post-processing
(e.g. constructing a parabolic surface for hertzian contact, computing the total
contact area). It is the user's responsibility to make the changes necessary,
but Tamaas provides some convenience functions to make this task easier. In the
module :py:mod:`mpi <tamaas._tamaas.mpi>`, the function :py:func:`local_shape
<tamaas._tamaas.mpi.local_shape>` gives the 2D shape of the local data if given
the global 2D shape (its counterpart :py:func:`global_shape
<tamaas._tamaas.mpi.global_shape>` does the exact opposite), while
:py:func:`local_offset <tamaas._tamaas.mpi.local_offset>` gives the offset of
the local data in the global :math:`x` dimension. These two functions mirror
FFTW's own data distribution `functions
<http://fftw.org/fftw3_doc/Basic-and-advanced-distribution-interfaces.html#Basic-and-advanced-distribution-interfaces>`_.
.. figure:: figures/mpi_data_distribution.svg
:align: center
:width: 75%
2D Data distribution scheme from FFTW. ``N0`` and ``N1`` are the
number of points in the :math:`x` and :math:`y` directions
respectively. The array ``local_N0``, indexed by the process rank,
give the local size of the :math:`x` dimension. The
:py:func:`local_offset <tamaas._tamaas.mpi.local_offset>` function
gives the offset in :math:`x` for each process rank.
The :py:mod:`mpi <tamaas._tamaas.mpi>` module also contains a function
:py:func:`sequential <tamaas._tamaas.mpi.sequential>` whose return value is
meant to be used as a context manager. Within the sequential context the default
communicator is ``MPI_COMM_SELF`` instead of ``MPI_COMM_WORLD``.
In addition to the low-level functions in :py:mod:`mpi <tamaas._tamaas.mpi>`,
many functions in Tamaas do not require special treatment, e.g. the statistics
functions.
Integration algorithm
---------------------
In its implementation of the volume integral operators necessary for
elastic-plastic solutions, Tamaas differenciates two way of computing the
intermediate integral along :math:`z` in the partial Fourier domain:
- Cutoff integration: because the intermediate integral involves kernels of the
form :math:`\exp(q(x-y))`, it is easy to truncate the integral when :math:`x`
and :math:`y` are far apart, especially for large values of :math:`q`. This
changes the complexity of the intermediate integral from
:math:`O(N_1N_2N_3^2)` (the naive implementation) to :math:`O(\sqrt{N_1^2 +
N_2^2}N_3^2)`.
- Linear integration: this method relies on a separation of variables
:math:`\exp(q(x-y)) = \exp(qx)\cdot\exp(-qy)`. This allows to break the
dependency in :math:`N_3^2` of the number of operations, so that the overall
complexity of the intermediate integral is :math:`O(N_1N_2N_3)`.
Details on both algorithms can be found in [1]_. Tamaas uses linear integration
by default because it is faster in many cases without introducing a truncation
error. Unfortunatly, it has a severe drawback when considering systems with a
fine surface discretization: due to :math:`q` increasing with the number of
points on the surface, the separated terms :math:`\exp(qx)` and
:math:`\exp(-qy)` may overflow and underflow respectively. Tamaas will warn
if that is the case, and users have two options to remedy the situation:
- Change the integration method by calling :func:`setIntegrationMethod
<tamaas._tamaas.Residual.setIntegrationMethod>` with the desired
:class:`integration_method <tamaas._tamaas.integration_method>` on the
:class:`Residual <tamaas._tamaas.Residual>` object you use in the computation.
- Compile Tamaas with the option ``real_type='long double'``. To make
manipulation of numpy arrays easier, a :class:`dtype` is provided in the
:py:mod:`tamaas` module which can be used to create numpy arrays compatible
with Tamaas' floating point type (e.g. ``x = np.linspace(0, 1,
dtype=tamaas.dtype)``)
Both these options negatively affect the performance, and it is up to the user
to select the optimal solution for their particular use case.
Computational methods
---------------------
Tamaas uses specialized numerical methods to efficiently solve elastic and
elastoplastic periodic contact problems. Using a boundary integral formulation
and a half-space geometry for the former allow (a) the focus of computational
power to the contact interface since the bulk response can be represented
exactly, (b) the use of the fast-Fourier transform for the computation of
convolution integrals. In conjunction with a boundary integral formulation of
the bulk state equations, a conjugate gradient approach is used to solve the
contact problem.
.. note:: The above methods are state-of-the-art in the domain of rough surface
contact. Below are selected publications detailing the methods used in
elastic contact with and without adhesion:
- Boundary integral formulation:
- Stanley and Kato (`J. of Tribology, 1997
<https://doi.org/10.1115/1.2833523>`_)
- Conjugate Gradient:
- Polonsky and Keer (`Wear, 1999
<https://doi.org/10.1016/S0043-1648(99)00113-1>`_)
- Rey, Anciaux and Molinari (`Computational Mechanics, 2017
<https://doi.org/10.1007/s00466-017-1392-5>`_)
- Frictional contact:
- Condat (`J. of Optimization Theory and Applications, 2012
<https://doi.org/10.1007/s10957-012-0245-9>`_)
For elastic-plastic contact, Tamaas uses a similar approach by implementing a
*volume* integral formulation of the bulk equilibrium equations. Thanks to
kernel expressions that are directly formulated in the Fourier domain, the
method reduces the algorithmic complexity, memory requirements and sampling
errors compared to traditional volume integral methods (Frérot, Bonnet, Anciaux
and Molinari, `Computer Methods in Applied Mechanics and Engineering, 2019
<https://doi.org/10.1016/j.cma.2019.04.006>`_, `arXiv:1811.11558
<https://arxiv.org/pdf/1811.11558.pdf>`_). The figure below shows a comparison
of run times for an elasticity problem (only a single solve step) between Tamaas
and `Akantu <https://gitlab.com/akantu/akantu>`_, a high-performance FEM code
using the direct solver `MUMPS <http://mumps.enseeiht.fr/>`_.
.. figure:: figures/complexity.svg
:align: center
:width: 75%
Comparison of run times between the volume integral implementation
(with cutoff integration) of Tamaas and an FEM solve step with a
Cholesky factorization performed by Akantu+MUMPS. :math:`N` is the
total number of points.
Further discussion about the elastic-plastic solver implemented in Tamaas can be
found in Frérot, Bonnet, Anciaux and Molinari, (`Computer Methods in Applied
Mechanics and Engineering, 2019 <https://doi.org/10.1016/j.cma.2019.04.006>`_,
`arXiv:1811.11558 <https://arxiv.org/pdf/1811.11558.pdf>`_).
.. [1] L. Frérot, “Bridging scales in wear modeling with volume integral
methods for elastic-plastic contact,” École Polytechnique Fédérale de
Lausanne, 2020 (Section 2.3.2). `10.5075/epfl-thesis-7640
<https://doi.org/10.5075/epfl-thesis-7640>`_
c4science · Help