Skip to content

Commit

Permalink
Merge branch 'master' into cgdna
Browse files Browse the repository at this point in the history
Conflicts:
	src/core/energy_inline.hpp
	src/core/interaction_data.cpp
	src/core/interaction_data.hpp
	src/features.def
	src/tcl/interaction_data_tcl.cpp
  • Loading branch information
fweik committed Dec 16, 2014
2 parents 3d32e2f + b86915b commit 03825e1
Show file tree
Hide file tree
Showing 177 changed files with 18,399 additions and 2,471 deletions.
580 changes: 495 additions & 85 deletions .cproject

Large diffs are not rendered by default.

7 changes: 7 additions & 0 deletions .gitignore
Expand Up @@ -7,9 +7,13 @@ config.status
src/acconfig.h.in
src/acconfig.h
aclocal.m4
*.m4
autom4te.cache
.deps/
myconfig.h
build/

*.*~

# created by configure
src/myconfig-final.h
Expand All @@ -18,6 +22,7 @@ doc/latexit.sh
myconfig-sample.h
testsuite/runtest.sh
tools/es_mpiexec
config/ltmain.sh

# object files and executable
*.o
Expand All @@ -26,6 +31,8 @@ src/libEspresso.a
Espresso
Espresso.install
build.eclipse
src/featuredefs.pyc
src/acconfig.hpp.in

# Espresso output
*.end
Expand Down
6 changes: 6 additions & 0 deletions .project
Expand Up @@ -5,6 +5,11 @@
<projects>
</projects>
<buildSpec>
<buildCommand>
<name>org.eclipse.cdt.autotools.core.genmakebuilderV2</name>
<arguments>
</arguments>
</buildCommand>
<buildCommand>
<name>org.eclipse.cdt.managedbuilder.core.genmakebuilder</name>
<triggers>clean,full,incremental,</triggers>
Expand Down Expand Up @@ -40,6 +45,7 @@
<nature>org.eclipse.cdt.managedbuilder.core.managedBuildNature</nature>
<nature>org.eclipse.cdt.managedbuilder.core.ScannerConfigNature</nature>
<nature>org.python.pydev.pythonNature</nature>
<nature>org.eclipse.cdt.autotools.core.autotoolsNatureV2</nature>
</natures>
<filteredResources>
<filter>
Expand Down
2 changes: 2 additions & 0 deletions NEWS
Expand Up @@ -11,6 +11,8 @@ New user-visible features
* New quartic bonded interaction
* New bonded coulomb interaction
* Lees Edwards boundaries for sheer stuff
* Espresso can now read pdb and gromacs topology files.
* Checkpointing for the correlator

Known bugs
----------
Expand Down
29 changes: 14 additions & 15 deletions configure.ac
Expand Up @@ -123,26 +123,15 @@ AX_CXX_VAR_PRETTYFUNC()
# set the optimization flags
# never overwrite users CXXFLAGS
if test "${CXXFLAGS+set}" != set; then
# test for -O5
AC_MSG_CHECKING([whether the compiler accepts -O5])
# test for -O3
AC_MSG_CHECKING([whether the compiler accepts -O3])
saved_CXXFLAGS=$CXXFLAGS
CXXFLAGS="-O5 $CXXFLAGS"
CXXFLAGS="-O3 $CXXFLAGS"
AC_COMPILE_IFELSE([AC_LANG_PROGRAM([],[])],[
AC_MSG_RESULT(yes);try_add_flag_res=yes
],[
AC_MSG_RESULT(no); CXXFLAGS=$saved_CXXFLAGS; try_add_flag_res=no
])
# if not, try O3 (e.g. llvm does not support O5)
AS_IF([test x$try_add_flag_res = xno], [
AC_MSG_CHECKING([whether the compiler accepts -O3])
saved_CXXFLAGS=$CXXFLAGS
CXXFLAGS="-O3 $CXXFLAGS"
AC_COMPILE_IFELSE([AC_LANG_PROGRAM([],[])],[
AC_MSG_RESULT(yes);try_add_flag_res=yes
],[
AC_MSG_RESULT(no); CXXFLAGS=$saved_CXXFLAGS; try_add_flag_res=no
])
],[])

##################################
# test for -Wall
Expand Down Expand Up @@ -375,6 +364,16 @@ AS_IF([test x$with_cuda != xno],[
LIBS=$save_noncuda_LIBS
AS_IF([test x$cuda_ok == xyes], [
# NVCC doesn't work with llvm's libc++ so we need to also
# compile the rest with libstdc++, if nvcc-compiled code is
# used
AC_MSG_CHECKING([whether the compiler accepts -stdlib=libstdc++])
saved_CXXFLAGS=$CXXFLAGS
CXXFLAGS="-stdlib=libstdc++ $CXXFLAGS"
AC_COMPILE_IFELSE([AC_LANG_PROGRAM([],[])],[ AC_MSG_RESULT(yes) ],[
AC_MSG_RESULT(no); CXXFLAGS=$saved_CXXFLAGS
])
AC_DEFINE(CUDA,[],[Whether CUDA is available])
],[
# reset standard compiler options
Expand Down Expand Up @@ -539,7 +538,7 @@ if test .$with_tk != .no; then
fi
# now test whether Tk can be found
if test .$with_tk = .yes; then
for version in $TK_VERSION tk8.5 tk8.4 tk8.3 tk8.2 tk; do
for version in $TK_VERSION tk8.6 tk8.5 tk8.4 tk8.3 tk8.2 tk; do
ES_ADDPATH_CHECK_LIB($version, Tk_Init, [use_tk=$version], [])
if test .$use_tk != .; then break; fi
done
Expand Down
49 changes: 47 additions & 2 deletions doc/ug/analysis-core.tex
Expand Up @@ -194,7 +194,9 @@ \subsubsection{Available observables}
\item \lit{density_profile} \var{particle\_specifications} \var{profile\_specifications} \\
Compute the density profile within the specified cube.
For profile specifications, see section~\ref{sec:DensProfSpec}.

\item \lit{force_density_profile} \var{particle\_specifications} \var{profile\_specifications} \\
Compute the force density profile within the specified cube.
For profile specifications, see section~\ref{sec:DensProfSpec}.
\item \lit{lb_velocity_profile} \var{particle\_specifications} \var{profile\_specifications} \\
Compute the Lattice-Boltzmann velocity profile within the specified cube.
For profile specifications, see section~\ref{sec:DensProfSpec}.
Expand All @@ -211,7 +213,11 @@ \subsubsection{Available observables}
\item \lit{lb_radial_velocity_profile} \\
Compute the Lattice-Boltzmann velocity profile in cylindrical coordinates.
For profile specifications, see section~\ref{sec:DensProfSpec}.
\end{itemize}
\item \lit{rdf} \var{type\_list type\_list [r\_min [r\_max [n\_bins]]]}\\
Compute the radial distribution function, see \ref{analyze:rdf}.
\item \lit{structure_factor} \var{order}\\
Compute the structure factor. Remember it scales as order$^3$, see \ref{analyze:structurefactor}.
\end{itemize}
The tclcommand observable is a helpful tool, that allows to make the
analysis framework much more versatile, by allowing
the evaluation of arbitrary tcl commands.
Expand Down Expand Up @@ -654,3 +660,42 @@ \subsection{The correlation algorithm: multiple tau correlator}
of the case of mean square displacement the difference is always positive,
resulting in a non-negligible systematic error. A more general
discussion is presented in Ref.~\cite{ramirez10a}.


\subsection{Checkpointing the correlator}
It is possible to checkpoint the correlator. Therby the data is written
directly to a file. It is probably usefull to write to a binary file, as this
is more precise.
The accuracy of the ascii files is not high.
\begin{essyntax}
\variant{1} correlation \var{id} write_checkpoint_binary \var{filename}
\variant{2} correlation \var{id} write_checkpoint_ascii \var{filename}
\end{essyntax}

In order to load a checkpoint, the correlator has to be initialized. Therefore
the observable(s) have to be created. Make sure that the correlator is exactly
initilized as it was when the checkpoint was created. If this is not
fullfilled, and e.g. the size of an observable has changed, loading the
checkpoint failes.
\begin{essyntax}
\variant{1} correlation \var{id} read_checkpoint_binary \var{filename}
\variant{2} correlation \var{id} read_checkpoint_ascii \var{filename}
\end{essyntax}
Depending on whether the checkpoint was written as binary or as text, the
corresponding variant for reading the checkpoint has to be used.

An simple example for checkpointing:
\begin{tclcode}
set pp [observable new particle_positions all]
set cor1 [correlation new obs1 $pp corr_operation square_distance_componentwise \
dt 0.01 tau_max 1000 tau_lin 16]
integrate 1000
correlation $cor1 write_checkpoint_binary "cor1.bin"
\end{tclcode}
And then to continue the simulation:
\begin{tclcode}
set pp [observable new particle_positions all]
set cor1 [correlation new obs1 $pp corr_operation square_distance_componentwise \
dt 0.01 tau_max 1000 tau_lin 16]
correlation $cor1 read_checkpoint_binary "cor1.bin"
\end{tclcode}
34 changes: 34 additions & 0 deletions doc/ug/citations.bib
Expand Up @@ -26,6 +26,40 @@ @INCOLLECTION{espresso2
address = {Berlin, Germany}
}

@article{Peskin2002,
author = {Peskin, Charles S},
title = {{The immersed boundary method}},
journal = {Acta Numerica},
year = {2003},
volume = {11},
pages = {479--517}
}

@article{Gompper1996,
author = {Gompper, G and Kroll, D M},
title = {{Random surface discretizations and the renormalization of the bending rigidity}},
journal = {Journal de Physique I},
year = {1996},
volume = {6},
pages = {1305--1320}
}

@article{Crowl2010,
author = {Crowl, Lindsay M and Fogelson, Aaron L},
title = {{Computational model of whole blood exhibiting lateral platelet motion induced by red blood cells}},
journal = {Int. J. Numer. Meth. Biomed. Engng.},
year = {2010},
volume = {26},
number = {3-4},
pages = {471--487},
}

@PHDTHESIS{KruegerThesis,
author = {Kr{\"u}ger, Timm},
title = {{Computer simulation study of collective phenomena in dense suspensions of red blood cells under shear}},
school = {Universit{\"a}t Bochum},
year = {2011}
}

@ARTICLE{elc,
author = {Axel Arnold and Jason {de Joannis} and Christian Holm},
Expand Down
76 changes: 76 additions & 0 deletions doc/ug/ibm.tex
@@ -0,0 +1,76 @@
% Copyright (C) 2010,2011,2012,2013,2014 The ESPResSo project
%
% This file is part of ESPResSo.
%
% ESPResSo is free software: you can redistribute it and/or modify it
% under the terms of the GNU General Public License as published by the
% Free Software Foundation, either version 3 of the License, or (at your
% option) any later version.
%
% ESPResSo 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
% General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see <http://www.gnu.org/licenses/>.
%

\chapter{Immersed Boundary Method for soft elastic objects}
\label{sec:ibm}

\begin{citebox}
Please contact the Biofluid Simulation and Modeling Group, University of Bayreuth, Germany if you plan to use this feature.
\end{citebox}

This section describes an alternative way to include soft elastic objects somewhat different from the previous chapter. In the Immersed Boundary Method (IBM), soft particles are considered as an infinitely thin shell filled with liquid (see e.g. ~\cite{Peskin2002, Crowl2010, KruegerThesis}). When the shell is deformed by an external flow it responds by elastic restoring forces which are transmitted into the fluid. In the present case, the inner and outer liquid are of the same type and are simulated using Lattice-Boltzmann.

Numerically, the shell is discretized by a set of marker points connected by triangles. The marker points are advected with \emph{exactly} the local fluid velocity, i.e., they do not possess a mass nor a friction coefficient (this is different from the Object-in-Fluid method of the previous chapter). We implement these marker points as virtual particles in \es which are not integrated using the usual velocity-verlet scheme, but instead are propagated using a simple Euler algorithm with the local fluid velocity (if the \texttt{IMMERSED_BOUNDARY} feature is turned on).

To compute the elastic forces, three new bonded interactions are defined ibm\_triel, ibm\_tribend and ibm\_volCons:

\begin{itemize}
\item ibm\_triel is a discretized elastic force with the following syntax
\begin{essyntax}
inter \var{ID} ibm_triel \var{ind1} \var{ind2} \var{ind3} \var{max} \var{law}
\begin{features}
\required{IMMERSED_BOUNDARY}
\end{features}
\end{essyntax}
where \var{ind1}, \var{ind2} and \var{ind3} represent the indices of the three marker points making up the triangle. The parameter \var{max}
specifies the maximum stretch above which the bond is considered broken. The final parameter \var{law} can be either
\begin{verbatim}
NeoHookean <k>
\end{verbatim}
or
\begin{verbatim}
Skalak <k1> <k2>
\end{verbatim}
which specifies the elastic law and its corresponding parameters (see e.g. \cite{KruegerThesis}).

\item ibm\_tribend is a discretized bending potential with the following syntax
\begin{essyntax}
inter \var{ID} ibm_tribend \var{ind1} \var{ind2} \var{ind3} \var{ind4} \var{method} \var{kb} \opt{flat|initial}
\begin{features}
\required{IMMERSED_BOUNDARY}
\end{features}
\end{essyntax}
where \var{ind1}, \var{ind2}, \var{ind3} and \var{ind4} are four marker points corresponding to two neighboring triangles. The indices \var{ind1} and \var{ind3} contain the shared edge. Note that the marker points within a triangle must be labelled such that the normal vector $\vec{n} = (\vec{r}_\text{ind2} - \vec{r}_\text{ind1}) \times (\vec{r}_\text{ind3} - \vec{r}_\text{ind1})$ points outward of the elastic object.

The parameter \var{method} allows to specify different numerical ways of computing the bending interaction. Currently, two methods are implemented, where the first one (\opt{Krueger}) follows \cite{KruegerThesis} and the second one (\opt{Gompper}) follows \cite{Gompper1996}. In both cases, \var{kb} is the bending modulus. The options \opt{flat} or \opt{initial} specify whether the reference shape is a flat configuration or whether the initial configuration is taken as reference shape, this option is only available for the \opt{Krueger} method.

\item ibm\_volCons is a volume-conservation force. Without this correction, the volume of the soft object tends to shrink over time due to numerical inaccuracies. Therefore, this implements an artificial force intended to keep the volume constant. If volume conservation is to be used for a given soft particle, the interaction must be added to every marker point belonging to that object. The syntax is
\begin{essyntax}
inter \var{ID} ibm_volCons \var{softID} \var{kv}
\begin{features}
\required{IMMERSED_BOUNDARY}
\end{features}
\end{essyntax}
where \var{softID} identifies the soft particle and \var{kv} is a volumetric spring constant \cite{KruegerThesis}.
\end{itemize}


%%% Local Variables:
%%% mode: latex
%%% TeX-master: "ibm"
%%% End:
57 changes: 56 additions & 1 deletion doc/ug/inter.tex
Expand Up @@ -1491,6 +1491,61 @@ \subsubsection{Additional P3M parameters}
\end{description}


\subsection{Coulomb Ewald GPU}
\label{sec:coulombewald}
\index{EwaldGPU method|mainindex}
\index{interactions!EwaldGPU|mainindex}

\begin{essyntax}
inter coulomb \var{l_B} ewaldgpu
\var{r_\mathrm{cut}} \alt{\var{K_\mathrm{cut}} \asep \{\var{K_\mathrm{cut,x}} \var{K_\mathrm{cut,y}} \var{K_\mathrm{cut,x}}\}} \var{alpha}
\begin{features}
\required{ELECTROSTATICS}
\end{features}
\end{essyntax}

This command activates the Ewald method to compute the electrostatic
interactions between charged particles. The far field is computed by the GPU with single precision and the near field by the CPU with double precision. It only works for the case of cubic boxes.
\begin{description}
\item[\var{l_B}] Bjerrum length as positive floating point number
\item[\var{r_\mathrm{cut}}] Real space cutoff as positive floating point number
\item[\var{K_\mathrm{cut}}] Reciprocal space cutoff as single positive integer
\item[\var{K_\mathrm{cut,xyz}}] Reciprocal space cutoff in x, y and z direction (relevant for noncubic boxes)
\item[\var{alpha}] Ewald parameter as positive floating point number
\end{description}

\subsubsection{Tuning Ewald GPU}
\label{ssec:tuneewaldgpu}
\begin{essyntax}
inter coulomb \var{l_B} ewaldgpu tune
accuracy \var{accuracy} precision \var{precision}
K_max \var{K_\mathrm{max}}
\begin{features}
\required{ELECTROSTATICS}
\end{features}
\end{essyntax}

The tuning algorithm first computes the optimal \var{r_\mathrm{cut}} and \var{alpha} for every \var{K_\mathrm{cut}} between one and \var{K_\mathrm{max}} as described in \cite{kolafa92}. Then the performance for all those (\var{K_\mathrm{cut}}, \var{r_\mathrm{cut}}, \var{alpha})-triplets will be measured via a short test simulation and the fastest will be chosen.

\begin{description}
\item[\var{accuracy}] Maximal allowed root mean square error regarding the forces
\item[\var{precision}] Determines how precise alpha will be computed
\item[\var{K_\mathrm{max}}] Maximal reciprocal space cutoff \var{K_\mathrm{cut}} to be tested in the tuning algorithm
\end{description}

\subsubsection{Tuning Alpha Ewald GPU}
\label{ssec:tunealphaewaldgpu}
\begin{essyntax}
inter coulomb \var{l_B} ewaldgpu tunealpha
\var{r_\mathrm{cut}} \alt{\var{K_\mathrm{cut}} \asep \{\var{K_\mathrm{cut,x}} \var{K_\mathrm{cut,y}} \var{K_\mathrm{cut,x}}\}} \var{precision}
\begin{features}
\required{ELECTROSTATICS}
\end{features}
\end{essyntax}

If \var{K_\mathrm{cut}} and \var{r_\mathrm{cut}} are given by the user, then \keyword{tunealpha} computes the optimal \var{alpha} with the chosen \var{precision} as described in \cite{kolafa92}. But in general \keyword{tune} should be chosen for tuning.


\subsection{Debye-H\"uckel potential}
\index{Debye-H\"uckel potential|mainindex}
\index{interactions!Debye-H\"uckel|mainindex}
Expand Down Expand Up @@ -2175,7 +2230,7 @@ \subsection{Tunable-slip boundary interaction}\label{sec:tunableSlip}
This method was tested for Dissipative Particle Dynamics but it is
intended for mesoscopic simulation methods in general. Note, that to
use tunable-slip boundary interactions you have to apply \textbf{two}
plane cell constraints with Lennard-Jones in addition to the
wall constraints with Lennard-Jones in addition to the
tunable-slip interaction. Make sure that the cutoff radius
\var{r_\mathrm{cut}} is larger than the cutoff radius of the
constraint Lennard-Jones interactions. Otherwise there is no
Expand Down

0 comments on commit 03825e1

Please sign in to comment.