Commit cf24164d authored by COULAUD Olivier's avatar COULAUD Olivier

Remoce Old Docs and add new FMath.hpp

parent 1ed80141
@article{fong09a,
author = {Fong, W and Darve, E},
doi = {DOI: 10.1016/j.jcp.2009.08.031},
file = {:Users/coulaud/Recherche/Bibliographie/Mendeley/Fong, Darve - 2009 - The black-box fast multipole method.pdf:pdf},
issn = {0021-9991},
journal = {Journal of Computational Physics},
keywords = {Chebyshev polynomials,Fast multipole method,Interpolation,Singular value decomposition},
mendeley-groups = {FMM},
number = {23},
pages = {8712--8725},
title = {{The black-box fast multipole method}},
url = {http://www.sciencedirect.com/science/article/pii/S0021999109004665},
volume = {228},
year = {2009}
}
@inproceedings{Warren1993,
address = {New York, NY, USA},
author = {Warren, M S and Salmon, J K},
booktitle = {Proceedings of the 1993 ACM/IEEE conference on Supercomputing},
doi = {http://doi.acm.org/10.1145/169627.169640},
file = {:Users/coulaud/Recherche/Bibliographie/Mendeley/Warren, Salmon - 1993 - A hashed Oct-Tree N-body algorithm.pdf:pdf},
isbn = {0-8186-4340-4},
mendeley-groups = {FMM},
pages = {12--21},
publisher = {ACM},
series = {Supercomputing '93},
title = {{A hashed Oct-Tree N-body algorithm}},
url = {http://doi.acm.org/10.1145/169627.169640},
year = {1993}
}
@Book{mason2003chebyshev,
author = {J C Mason and D C Handscomb},
title = {Chebyshev polynomials},
publisher = {Chapman \& Hall/CRC},
year = {2003}
}
This diff is collapsed.
\section{Splitting the outer product in the interpolation polynomial}
\label{sec:splitouter}
We introduce the notation of the interpolation polynomial based on the
notation from \cite{mason2003chebyshev}, the interpolation of a function
$p(x)$ reads as
\begin{equation}
\label{eq:masondef}
p(x) \sim \sum_{n=0}^{\ell-1} \gamma_n c_n T_n(x)
\end{equation}
with $\gamma_0=1$ and $\gamma_n=2$ for $n\ge1$ and with
\begin{equation}
\label{eq:chebcoeff}
c_n = \frac{1}{\ell} \sum_{k=1}^\ell p(\bar x_n) T_k(\bar x_n)
\end{equation}
After inserting Eqn.~(\ref{eq:chebcoeff}) into Eqn.~(\ref{eq:masondef}) we end
up with
\begin{equation}
\label{eq:outprod}
p(x) \sim \frac{1}{\ell} \sum_{n=1}^\ell p(\bar x_n) \sum_{k=0}^{\ell-1}
\gamma_k T_k(\bar x_n) T_k(x)
\end{equation}
In the discrete case, i.e., we have $N$ scattered points $\{x_i\}_{i=1}^N
\subset [-1,1]$, we write Eqn.~\eqref{eq:outprod} in matrix notation as
\begin{equation}
\label{eq:interpolmatnot}
\Mat{p} \sim \Mat{T}_x \Mat{\bar T}_x^\top \Mat{\bar p}
\end{equation}
with the matrices $\Mat{T}_x \in \mathbb{R}^{N\times\ell}$ with
$(\Mat{T}_x)_{ik} = T_{k-1}(x_i)$ and $\Mat{\bar T}_x \in
\mathbb{R}^{\ell\times\ell}$ with $(\Mat{\bar T}_x)_{nk} =
\nicefrac{\gamma_{k-1}}{\ell} \, T_{k-1}(\bar x_n)$. Using this notation the
interpolated kernel function $K(x,y)$ reads as
\begin{equation}
\label{eq:kinterpolmatnot}
\Mat{K} \sim \Mat{T}_x \Mat{\bar T}_x^\top \Mat{\bar K} (\Mat{T}_y
\Mat{\bar T}_y^\top)^\top
\end{equation}
With the low-rank representation $\Mat{\bar K} \sim \Mat{UV}^T$ we compute
$\Mat{\bar U} = \Mat{\bar T}_x^\top \Mat{U}$ and $\Mat{\bar V} = \Mat{\bar
T}_y^\top \Mat{V}$ and write Eqn.~\eqref{eq:kinterpolmatnot} as
\begin{equation}
\label{eq:kinterpolmatnot1}
\Mat{K} \sim \Mat{T}_x \Mat{\bar U} \Mat{\bar V}^\top \Mat{T}_y^\top
\end{equation}
\paragraph{Can we still use symmetries?} Recall the approach we proposed in
Sec.~\ref{sec:m2l}. It exploits the symmetries in the arrangement of the
far-field interactions. As we see from Eqn.~\eqref{eq:permut} in that case we
need permutation matrices $\Mat{P}_t$ for the $t$-th far-field interaction and
Eqn.~\eqref{eq:kinterpolmatnot} becomes
\begin{equation}
\label{eq:symsplit}
\Mat{K}_t \sim \Mat{T}_x \Mat{\bar T}_x^\top (\Mat{P}_t \Mat{\bar K}_{p(t)}
\Mat{P}_t^\top) (\Mat{T}_y \Mat{\bar T}_y^\top)^\top.
\end{equation}
With $\Mat{\bar K}_{p(t)} \sim \Mat{U}_{p(t)} \Mat{V}_{p(t)}^\top$ the
matrices in the outer product in Eqn.~\eqref{eq:kinterpolmatnot1} become
$\Mat{\bar U}_t = \Mat{\bar T}_x^\top \Mat{P}_t \Mat{U}_{p(t)}$ and $\Mat{\bar
V}_t = \Mat{\bar T}_y^\top \Mat{P}_t \Mat{V}_{p(t)}$. Thus, it is not
possible anymore to use the fact that due to symmetries we can express all
$316$ far-field interactions by permutations of $16$ only.
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "main"
%%% End:
This diff is collapsed.
@article{Schmidt1991,
abstract = {The Rokhlin-Greengard fast multipole algorithm for evaluating Coulomb and multipole potentials has been implemented and analyzed in three dimensions. The implementation is presented for bounded charged systems and systems with periodic boundary conditions. The results include timings and error characterizations},
author = {Schmidt, K. E. and Lee, Michael a.},
doi = {10.1007/BF01030008},
file = {:Users/coulaud/Recherche/Bibliographie/Mendeley/Schmidt, Lee - 1991 - Implementing the fast multipole method in three dimensions.pdf:pdf},
issn = {0022-4715},
journal = {Journal of Statistical Physics},
keywords = {Ewald summation,Periodic boundary conditions,fast multipole method,macro summation,many-body problem,n-body,spherical harmonique},
mendeley-groups = {FMM},
month = jun,
number = {5-6},
pages = {1223--1235},
title = {{Implementing the fast multipole method in three dimensions}},
url = {http://link.springer.com/10.1007/BF01030008},
volume = {63},
year = {1991}
}
@article{Heyes1981,
abstract = {The Coulomb potentials and fields in infinite point charge lattices are represented by series expansions in real and reciprocal space using a method similar to that devised by Bertaut. A systematic investigation is made into the relative convergence characteristics of series derived by varying the charge spreading function which is required in the theory. A number of formulas which are already in use are compared with new expressions derived by this method.},
author = {Heyes, D. M.},
doi = {10.1063/1.441285},
file = {:Users/coulaud/Recherche/Bibliographie/Mendeley/Heyes - 1981 - Electrostatic potentials and fields in infinite point charge lattices.pdf:pdf},
issn = {00219606},
journal = {The Journal of Chemical Physics},
keywords = {Electrostatics,Lattice theory},
mendeley-groups = {Chimie},
mendeley-tags = {Electrostatics,Lattice theory},
number = {3},
pages = {1924},
title = {{Electrostatic potentials and fields in infinite point charge lattices}},
url = {http://link.aip.org/link/?JCP/74/1924/1\&Agg=doi},
volume = {74},
year = {1981}
}
@article{DeLeeuw1980,
author = {{De Leeuw}, S W and Perram, J W and Smith, E R},
file = {:Users/coulaud/Recherche/Bibliographie/Mendeley/De Leeuw, Perram, Smith - 1980 - Simulation of electrostatic systems in periodic boundary conditions. I. Lattice sums and dielectric con.pdf:pdf},
journal = {Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences},
keywords = {Periodic Boundary Conditions,dieclectric constant,lattice sums},
mendeley-groups = {Chimie},
mendeley-tags = {Periodic Boundary Conditions,dieclectric constant,lattice sums},
number = {1752},
pages = {27--56},
title = {{Simulation of electrostatic systems in periodic boundary conditions. I. Lattice sums and dielectric constants}},
url = {http://www.jstor.org/stable/2397238 http://rspa.royalsocietypublishing.org/content/373/1752/27.short},
volume = {373},
year = {1980}
}
\documentclass[12pt,letterpaper,titlepage]{article}
\usepackage{algorithm2e}
\usepackage{listings}
\usepackage{geometry}
\usepackage{graphicx}
\usepackage{tikz}
\usepackage[hypertexnames=false, pdftex]{hyperref}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% use:$ pdflatex TreeBuilder.tex
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\author{Berenger Bramas, Olivier Coulaud, Cyrille Piacibello}
\title{ScalFMM - Tree Building (Draft)}
\date{\today}
%% Package config
\lstset{language=c++, frame=lines}
\RestyleAlgo{boxed}
\geometry{scale=0.8, nohead}
\hypersetup{ colorlinks = true, linkcolor = black, urlcolor = blue, citecolor = blue }
%% Remove introduction numbering
\setcounter{secnumdepth}{-1}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{document}
\maketitle{}
\subsection*{How to Use}
There is a ccp class, named FTreeBuilder, which can fill a tree with
particles.
\subsubsection*{Old way}
The old method is the following :
\begin{enumerate}
\item Load a Particle File.
\item Create an octree.
\item Insert each particle in the Tree.
\end{enumerate}
Code will looks like that :
\begin{lstlisting}[frame=single]
FFmaGenericLoader loader(filename);
OctreeClass tree(TreeHeight, SubTreeHeight,
loader.getBoxWidth(), loader.getCenterOfBox());
int nbParts = loader.getNumberOfParticles();
for(int idxPart = 0 ; idxParts < nbParts ; ++idxParts){
FPoint position;
FReal physicalValue;
loader.fillParticles(&position,&physicalValue);
tree.insert(position,physicalValue);
}
\end{lstlisting}
\subsubsection*{New way}
\begin{enumerate}
\item Load a Particle File.
\item Create an octree.
\item Store each part in a container.
\item Call TreeBuilder on that container and a pointer to the tree.
\end{enumerate}
\newpage
Code will looks like that :
\begin{lstlisting}[frame=single]
FFmaGenericLoader loader(filename);
OctreeClass tree(TreeHeight, SubTreeHeight,
loader.getBoxWidth(), loader.getCenterOfBox());
int nbParts = loader.getNumberOfParticles();
FmaRWParticle * parts = new FmaRWParticle[nbParts];
for(int idxPart = 0 ; idxParts < nbParts ; ++idxParts){
FPoint position;
FReal physicalValue;
loader.fillParticles(&position,&physicalValue);
parts[idxPart].setPosition(position);
parts[idxPart].setPhysicalValue(physicalValue);
}
typedef FTreeBuilder<FmaRWParticle,OctreeClass,LeafClass> TreeBuilder;
TreeBuilder::BuildTreeFromArray(parts,nbParts,&tree);
\end{lstlisting}
\subsection*{How it works}
The main objective is to reduce the time wasted in the insertion. The
operation is long for two reasons, the fact that the parts are
included one after the other, and the fact that we must go through all
the levels of the tree to find where to store the part.
Finally, the parts are not stored as an array of struct, but as
several arrays for each particle's attribute. Those array need to be
redimensionned each time we had a particle.
Main Idea to improve the time used : Insert arrays of particles
instead of particles.
The ParticleClass (in the example, it's a FmaRWParticle\textless
R,W\textgreater)used need to provide at least a getPosition() and a
getPhysicalValue() method.
Step by step analysis :
\begin{enumerate}
\item The array is copied into an array of IndexedParticle structure :
\begin{lstlisting}[frame=single]
struct IndexedParticle{
public:
MortonIndex index;
ParticleClass particle;
operator MortonIndex() const {
return this->index;
}
bool operator<=(const IndexedParticle& rhs){
return this->index <= rhs.index;
}
};
\end{lstlisting}
The index are setted using GetTreeCoordinate() (method copied from
OctreeClass) and GetMortonIndex() methods from FTreeCoordinate.
\item The array is then sorted in order to know where each particle
belongs. This step can be skipped if the array is already sorted
along MortonIndex. To skip it, just pass true as the third args of
BuildTreeFromArray();
If the array needs to be sorted, then it is, using FQuickSort.hpp,
our threaded implementation of the quick sort algorithm.
\item After the sort, we count the number of different leaves. And
we store the offset in an other array. Furthermore, we copy for
each leaves the position into a FPoint array, and the physical
value into a FReal array.
Example of MortonIndex Distributions and Offset array corresponding :
Parts :
\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|}
\hline
0&1&1&2&4&4&4&5&6&6\\
\hline
\end{tabular}
Offset :
\begin{tabular}{|c|c|c|c|c|c|}
\hline
0&1&3&4&7&8\\
\hline
\end{tabular}
\item We create every leaves in the tree (note that there are
empty). And we store the pointer to each leaves in an array.
\item For each leaves created, we fill it with a new method in
FSimpleParticle and FBasicParticleContainer, pushArray. This
method is called on the array of position and the array of
physical value.
\end{enumerate}
\end{document}
\ No newline at end of file
\documentclass{minimal}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{tikz}
\usetikzlibrary{trees}
\tikzstyle{coordaxis}=[draw=red!50!black, thick, ->]
\tikzstyle{mapping}=[draw=blue!50!black, thick, ->]
\tikzstyle{connect}=[draw=green!50!black, thick, ->]
\newcommand{\sR}{\mathbb{R}}
\begin{document}
% grow=down,
% every node/.style={draw, circle, thin},
% edge from parent/.style={-latex, thick, draw}
\begin{tikzpicture}[scale=0.6]
\draw[fill=blue!20] (0,8) rectangle +(8,8);
\draw[fill=blue!20] (8,0) rectangle +(8,8);
\draw[fill=blue!20] (0,0) rectangle +(4,4);
\draw[fill=blue!20] (12,12) rectangle +(4,4);
\draw[fill=brown!50] (4,0) rectangle +(4,4);
\draw[fill=brown!50] (0,4) rectangle +(4,4);
\draw[fill=brown!50] (4,4) rectangle +(4,4);
\draw[fill=brown!50] (8,8) rectangle +(4,4);
\draw[fill=brown!50] (8,12) rectangle +(4,4);
\draw[fill=brown!50] (12,8) rectangle +(4,4);
%Draw the initial octree
\draw[step = 8] (0,0) grid (16,16);
\draw[step = 4] (8,8) grid(16,16);
\draw[step = 4] (0,8) grid(8,16);
\draw[step = 2] (0,4) grid(4,8);
\draw[step = 2] (4,4) grid (8,8);
\draw[step = 2] (4,0) grid (8,4);
\draw[step = 1] (2,4) grid (4,8);
\draw[step = 1] (0,6) grid (2,8);
\draw[step = 0.5] (2,6) grid (3,7);
\draw[step = 2] (8,8) grid(12,16);
\draw[step = 2] (12,8) grid (16,12);
\draw[step = 1] (12,8) grid (14,12);
\draw[step = 1] (14,8) grid (16,10);
\draw[step = 0.5] (13,9) grid (14,10);
\end{tikzpicture}
%sibling distance=2cm,
\begin{tikzpicture}[level distance=1.5cm, grow=down,
every node/.style={draw, circle, thin},
edge from parent/.style={-latex, thick, draw},
level 1/.style={sibling distance=5cm},
level 2/.style={sibling distance=2cm},
level 3/.style={sibling distance=1.cm}
]
\node (P) {Root}
child {node [fill=blue] (C1) {}
child {node (T) {}
child { node (TTT) {}
child {node (TTTT) {}
child {node (TTTTT) {}
}
}
}
}
child {node (U) {}
child { node [fill=blue] (TTT) {}
child {node (TTTT) {}
child {node (TTTTT) {}
}
}
child {node [fill=blue] (TTTTU) {}
child {node (TTTTU1) {}}
child {node (TTTTU2) {}}
}
}
}
}
child {node (C2) {}}
% child {node (C3) {}}
child {node [fill=blue, text=white] (C4) {$j$}
child {node (T) {}
child {node (TT) {}
child { node [fill=blue] (TTT) {}
child {node (TTTT) {}}
child {node (TTTU) {}}
}
}
}
child {node (U) {}}
};
\path (P) -- coordinate[midway] (PQ) (C1);
\path (P) -- coordinate[midway] (PR) (C2);
%\draw (PQ) to[bend right=22] (PR);
\end{tikzpicture}
\begin{tikzpicture}
[level distance=30mm, %level/.style={sibling distance=36mm/#1}]
every node/.style={draw, circle, thin},
% edge from parent/.style={-latex, thick, draw},
level 1/.style={sibling distance=12cm},
level 2/.style={sibling distance=30mm},
level 3/.style={sibling distance=7mm},
scale=0.5, rotate=90 ]
\coordinate
child foreach \x in {0,1,2,3} {node{}
{child foreach \y in {0,1,2,3} {node{}
{child foreach \z in {0,1,2,3} {node{}}} }}};
\end{tikzpicture}
\begin{tikzpicture}
%\node {root} child [red] foreach \name in {1,2} {node {\name}}
\end{tikzpicture}
\end{document}
\ No newline at end of file
\documentclass[10pt]{article}
%
\usepackage{a4wide}
\usepackage{amsmath}
\usepackage{graphics,subfig}
\usepackage{tikz}
\usepackage{pgfkeys}
%
\usetikzlibrary{shapes}
% p 48
\title{Distribution of particles into ScalFMM}
\author{O. Coulaud}
\date{March, 29th 2014\\\vspace{1cm}Version 0.1}
\begin{document}
\maketitle
\tableofcontents
\newpage
\section{Uniform distribution for surface}
For parametric surface, J Williamson in \cite{Williamson} proposed a method to construct uniform distribution of points on such surface. For ellipsoid, this method was improved in \cite{ChenGlotzer}. Here we follow their approach to construct uniform distribution.
\subsection{Sphere distribution}
Consider that we want $N$ point uniformly set on the surface of the unit sphere, the natural choice is to sample uniformly the angles $\theta$ and $\phi$ of the spherical coordinates. Unfortunately, such choice lead to a concentration of points near the pole. As a surface area is given by $sin\phi\,d\phi\,d\theta = -d(cos\phi)\,d\theta$, we will sample $cos\phi$ rather than $\Phi$. So, we choose $u$ and $v$ two random variable on [0,1] and
\begin{eqnarray*}
\theta &=& 2\pi u,\\
\phi &=& cos^{-1}(2v-1).
\end{eqnarray*}
then we obtain a uniform distribution of points on the unit sphere.
\begin{figure}[ht]
\centering
\includegraphics[width=0.3\textwidth]{unitsphere}
\caption{$5\,000$ points distribution on unit sphere.}%
\end{figure}
\texttt{generateDistributions -unitsphere -N 5000 -fout unitsphere.fma -fvisuout unitsphere.vtp }
\subsection{Ellipsoid distribution}
Here we want to construct a uniform distribution on an ellipsoid defined by the equation
\begin{equation}
f(x,y,z) = (\frac{x}{a})^2+(\frac{y}{b})^2+(\frac{z}{c})^2-1,
\end{equation}
and the parameterization
\begin{eqnarray*}
x &=& a\cos(u)\cos(v),\\
y &=& b\cos(u)\sin(v), \\
z &=& c\sin(u).
\end{eqnarray*}
with $u\in [0, 2\pi[$ and $v\in [0, \pi[$. if you sample $u$ and $v$ then we obtain a non uniform distribution on the ellipsoid see section~ \ref{nonuniEllipsoid}.
\subsubsection{prolate}
A prolate spherical geometry is an ellipsoid where $a=b$ and $c>a$. On the case $g_{max}=1+a^2$ is obtain when $z=o$. Then to construct a uniform distribution we proceed as follow
\begin{description}
\item[step 1] Build a uniform point on the sphere and map it on the prolate surface. Let u a random point in $[-1,1]$ and v in $[0,\pi]$, then
\begin{eqnarray*}
x &=& a\sqrt{1-u^2}\cos(v),\\
y &=& b\sqrt{1-u^2}\sin(v), \\
z &=& c u.
\end{eqnarray*}
\item[step 2] Correct the distribution. Choose a random point $\xi$ in $[0,1]$ and if
$$ x^2+y^2 + \frac{a^4}{c^4} z^2 < a^2 \xi^2 $$
is true then we keep the point otherwise we reject it.
\end{description}
\begin{figure}[h]
\centering
\begin{minipage}{0.45\textwidth}%
\includegraphics[width=1.0\textwidth]{ellipsoidUnif}
\caption{$20\,000$ points distribution on 2:2:4 ellipsoid. Less than $15 \%$ of points has been rejected.}%
\end{minipage}%
\qquad
\begin{minipage}{0.45\textwidth}%
\includegraphics[width=1.2\textwidth]{prolate}
\caption{$20\,000$ points distribution on 1:1:10 ellipsoid and only $21\%$ of the tested points in step 2 has been rejected.}%
% \label{fig:2figsB}%
\end{minipage}%
\end{figure}
\texttt{generateDistributions -prolate -ar 1:1:10 -N 20000 -fout prolate.bfma -fvisuout prolate.vtp}
\section{Non uniform distribution}
\subsection{Ellipsoid}
\label{nonuniEllipsoid}
The simplest way to sample and ellipsoid of size a:b:c is to consider its spherical coordinates. Let $u$ and $v$ two random numbers between 0 and 1. After mapping u (resp. v) in $[-\pi/2,\pi/2]$ resp. ($[-\pi,\pi]$), the Cartesian coordinates are
\begin{eqnarray*}
x &=& a\cos(u)\cos(v),\\
y &=& b\cos(u)\sin(v), \\
z &=& c\sin(u).
\end{eqnarray*}
As shown in the figure~\ref{Fig-nonUnifEllipsoid} we obtain a concentration of points near the poles $(0,0,\pm c)$.
\begin{figure}[hbt]
\centering
\includegraphics[width=0.5\textwidth]{ellipsoid}
\caption{$20\,000$ points distribution on an ellipsoid of aspect ratio 2:2:4.}%
\label{Fig-nonUnifEllipsoid}
\end{figure}
This bias on the pole could be reduced by choosing the same approach that we use to build a uniform distribution on the unit sphere?
\texttt{generateDistributions -ellipsoid -ar 2:2:4 -N 20000 -fout ellipsoid.bfma -fvisuout ellipsoid.vtp}
If you consider the
\subsection{Plummer Model}
This is a hard test case in astrophysics problem, and it models a globular cluster of stars, which is highly non uniform. It is called the plummer distribution. To construct such distribution, first we construct a uniform points distribution on the unit sphere. Second, the radius is chosen according to the plummer distribution (double power law in astrophysics). We consider $u$ a random number between 0 and 1, then the associated radius is given by
\begin{equation*}
r = 1.0/\sqrt{u^{-2/3}-1},
\end{equation*}
and the total mass is one. Then, $m_i = \frac{1}{Npt}$.
\begin{figure}[h]
\centering
\begin{minipage}{0.45\textwidth}%
\includegraphics[scale=0.37,angle=-90]{plummerHistogramme}
\caption{Radius distribution}%
\end{minipage}%
\qquad
\begin{minipage}{0.45\textwidth}%
\includegraphics[width=1.0\textwidth]{plummer3D}
\caption{$50\,000$ point distribution.}%
\end{minipage}%
\end{figure}
The command to generate such distribution is\\
\texttt{generateDistributions -plummer -radius 10 -N 50000 -fout plummer.bfma -fvisuout plummer.vtp
}
The Plummer 3-dimensional density profile is
\begin{equation}
\rho_P(r) = \frac{3 M}{4\pi a^3} (1+\frac{r^2}{a^2})^{-\frac{5}{2}}
\end{equation}
where M is the total mass of the cluster and $a$ the Plummer radius.
The corresponding potential is
\begin{equation}
\Phi_P(r) = - \frac{G M}{\sqrt{r^2+a^2}}
\end{equation}
In N-body units, $G = M = 1$ and $a = 3\pi/16 \sim 0.589$
\subsection{Diagonal Model}
%, shape end size=.5cm},decoration={shape start size=.5cm, shape end size=.125cm
\tikzset{bigger/.style={decoration={shape start size=.125cm}}}
\tikzset{myblack/.style={draw,fill=black}}
\begin{tikzpicture}[scale=0.5]
Lmax = 2;
a := 10;
\foreach \l in {1,2,3,4,5} {
\draw (10-20/2^\l,10-10/2^\l) -- (10.0,10-10/2^\l) ;
\draw (10-10/2^\l,10-20/2^\l) -- (10-10/2^\l,10.0) ;
\draw [myblack] (10-15/2^\l,10-15/2^\l) circle (0.1) ;
}
\draw [myblack] (10-5/2^5,10-5/2^5) circle (0.1) ;
\draw (0,0) rectangle (10,10);
\end{tikzpicture}
% \node at (10-15/2^\l,10-15/2^\l) [shape=circle,size=6mm,draw=black,fill=black] {} ;
\begin{tikzpicture}[scale=0.5]
Lmax = 2;
a := 10;
\foreach \l in {1,2,3,4,5} {
\draw (10-20/2^\l,10-10/2^\l) -- (10.0,10-10/2^\l) ;
\draw (10-10/2^\l,10-20/2^\l) -- (10-10/2^\l,10.0) ;
\draw [myblack] (10-15/2^\l,10-15/2^\l) circle (0.1) ;
}
\draw [myblack] (10-5/2^5,10-5/2^5) circle (0.1) ;
\draw (0,0) rectangle (10,10);
\end{tikzpicture}
\newpage
\begin{thebibliography}{1}
\bibitem{Williamson} J. F. Williamson, “Random selection of points distributed on curved surfaces,” Physics in Medicine and Biology, vol. 32, no. 10, pp. 1311–1319, Oct. 1987.
\bibitem{ChenGlotzer} T. Chen and S. C. Glotzer, “Simulation studies of a phenomenological model for elongated virus capsid formation,” Physical review. E, Statistical, nonlinear, and soft matter physics, vol. 75, pp. 1–25, 2007.
\end{thebibliography}
\newpage
\section{Annexe}
\subsection{Gnuplot script for histogram}
\begin{verbatim}
clear
reset
set key off
set border 3
# Add a vertical dotted line at x=0 to show centre (mean) of distribution.
set yzeroaxis
# Each bar is half the (visual) width of its x-range.
set boxwidth 0.05 absolute
set style fill solid 1.0 noborder
bin_width = 0.1;
bin_number(x) = floor(x/bin_width)
rounded(x) = bin_width * ( bin_number(x) + 0.5 )
set xlabel "Radius"
set ylabel "Number of Particles"
#set terminal postscript enhanced color 'Helvetica' 20
#set output 'plummerHistogramme.eps'
plot 'plummerNewSort.txt' using (rounded($1)):(1) smooth frequency with boxes
\end{verbatim}
\end{document}
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
......@@ -328,6 +328,24 @@ struct FMath{
return 1.0/sqrt(inValue);
}
#ifdef SCALFMM_USE_SSE
static __m128 Exp(const __m128 inV){
float ptr[4];
_mm_storeu_ps(ptr, inV);
for(int idx = 0 ; idx < 4 ; ++idx){
ptr[idx] = std::exp(ptr[idx]);
}
return _mm_loadu_ps(ptr);
}
static __m128d Exp(const __m128d inV){