Nuclearistes/algorithms.org
2020-12-12 01:15:06 +01:00

486 lines
20 KiB
Org Mode

#+TITLE: Important algorithms for CIPSI
#+DATE: 14/12/2020
#+AUTHOR: Anthony Scemama
#+startup: beamer
#+LaTeX_HEADER: \institute{Laboratoire de Chimie et Physique Quantiques, IRSAMC, UPS/CNRS, Toulouse}
#+LATEX_CLASS: beamer
#+LaTeX_CLASS_OPTIONS:[aspectratio=169]
#+BEAMER_THEME: pterosor
#+LaTeX_HEADER: \usepackage{minted}
#+LaTeX_HEADER: \usepackage[utf8]{inputenc}
#+LaTeX_HEADER: \usepackage[T1]{fontenc}
#+LaTeX_HEADER: \usepackage{hyperref}
#+LaTeX_HEADER: \usepackage{mathtools}
#+LaTeX_HEADER: \usepackage{physics}
#+LaTeX_HEADER: \definecolor{darkgreen}{rgb}{0.,0.6,0.}
#+LaTeX_HEADER: \definecolor{darkblue}{rgb}{0.,0.2,0.7}
#+LaTeX_HEADER: \definecolor{darkred}{rgb}{0.6,0.1,0.1}
#+LaTeX_HEADER: \definecolor{darkpink}{rgb}{0.7,0.0,0.7}
#+EXPORT_EXCLUDE_TAGS: noexport
#+OPTIONS: H:2 num:t toc:nil \n:nil @:t ::t |:t ^:t -:t f:t *:t <:t
#+OPTIONS: TeX:t LaTeX:t skip:nil d:nil todo:t pri:nil tags:not-in-toc
#+LATEX: \newcommand{\mcenter}[1]{\multicolumn{1}{c}{#1}}
#+LATEX: \newcommand{\Ndet}{{\textcolor{darkgreen}{N_\text{det}}}}
#+LATEX: \newcommand{\Evar}{\textcolor{darkgreen}{E_\text{var}}}
#+LATEX: \newcommand{\Ept}{\textcolor{red}{E_\text{PT2}}}
#+LATEX: \newcommand{\de}{\textcolor{red}{\delta E(\alpha)}}
#+LATEX: \newcommand{\mH}{\mathcal{H}}
#+LATEX: \newcommand{\mDup}{{\textcolor{orange}{D_k^\uparrow}}}
#+LATEX: \newcommand{\mDupp}{{\textcolor{orange}{D_{k'}^\uparrow}}}
#+LATEX: \newcommand{\mDdn}{{\textcolor{orange}{D_m^\downarrow}}}
#+LATEX: \newcommand{\mDdnn}{{\textcolor{orange}{D_{m'}^\downarrow}}}
#+LATEX: \newcommand{\Ndetup}{{\textcolor{orange}{N_\text{det}^\uparrow}}}
#+LATEX: \newcommand{\Ndetdn}{{\textcolor{orange}{N_\text{det}^\downarrow}}}
#+LATEX: \newcommand{\mD}{{\textcolor{darkgreen}{\mathcal{D}}}}
#+LATEX: \newcommand{\mPsi}{{\textcolor{darkgreen}{\Psi}}}
#+LATEX: \newcommand{\mA}{\textcolor{red}{\mathcal{A}}}
#+LATEX: \newcommand{\ma}{{\textcolor{red}{\alpha}} }
#+LATEX: \newcommand{\mi}{{\textcolor{darkgreen}{I}} }
#+LATEX: \newcommand{\mpi}{{\textcolor{darkgreen}{p_I}} }
#+LATEX: \newcommand{\mci}{{\textcolor{darkgreen}{c_I}} }
#+LATEX: \newcommand{\mcj}{{\textcolor{darkgreen}{c_J}} }
#+LATEX: \newcommand{\epsi}{{\textcolor{darkgreen}{\epsilon}_\mi} }
#+LATEX: \newcommand{\mj}{{\textcolor{darkgreen}{J}} }
#+LATEX: \newcommand{\nea}{{N_\text{elec}^\uparrow}}
#+LATEX: \newcommand{\neb}{{N_\text{elec}^\downarrow}}
#+LATEX: \newcommand{\nmo}{{N_\text{MO}}}
#+LATEX: \newcommand{\Nsamples}{N_\text{samples}}
#+LATEX: \newcommand{\uniqueness}{{\textcolor{blue}{uniqueness}} }
#+LATEX: \newcommand{\mm}{{\textcolor{blue}{m}} }
#+LATEX: \newcommand{\mk}{{\textcolor{blue}{k}} }
#+LATEX: \newcommand{\mM}{{\textcolor{blue}{M}} }
#+LATEX: \newcommand{\epsik}{{\textcolor{darkgreen}{\epsilon_{I_{\textcolor{blue}{k}}}} }}
* Efficient direct CI
** Sorting
*** Popular misconception *Sorting is /not/ $\mathcal{O}(N \log(N))$* :
#+LATEX: \textcolor{darkgreen}{sorting is $\mathcal{O}(N \log(M))$}
(linear in $N$, log in $M$)
*** . :B_ignoreheading:
:PROPERTIES:
:BEAMER_env: ignoreheading
:END:
- A is an array of $N$ integer values
- The bitmask is an integer with only one bit set to one (00001000)
#+HEADER: :includes (list "<stddef.h>" "<stdio.h>")
#+BEGIN_SRC C
void radix_sort(int* A, size_t N, int bitmask) {
if (bitmask == 0) return;
int left[N], right[N];
int p=0 ; int q=0 ;
for (int i=0 ; i<N ; i++) {
if (A[i] & bitmask) { right[q] = A[i]; q++; }
else { left [p] = A[i]; p++; } }
radix_sort(left , p, bitmask >> 1 ) ;
radix_sort(right, q, bitmask >> 1 ) ;
for (int i=0 ; i<p ; i++) { A[ i] = left [i] ; }
for (int i=0 ; i<q ; i++) { A[p+i] = right[i] ; }
}
#+END_SRC
*** No export :noexport:
#+BEGIN_SRC C :export none
int main() {
int A[14] = {5,2,4,5,7,8,9,65,4,3,3,5,6,72};
for (int i=0 ; i<14 ; i++)
printf("%d ", A[i]);
printf("\n");
radix_sort(A,14, 1 << 24);
for (int i=0 ; i<14 ; i++)
printf("%d ", A[i]);
printf("\n");
}
#+END_SRC
#+RESULTS:
| 5 | 2 | 4 | 5 | 7 | 8 | 9 | 65 | 4 | 3 | 3 | 5 | 6 | 72 |
| 2 | 3 | 3 | 4 | 4 | 5 | 5 | 5 | 6 | 7 | 8 | 9 | 65 | 72 |
** Double-determinant
A Slater determinant can be written as a Waller-Hartree double determinant
\[ |I \rangle = \hat{I}\, |\rangle =
-1^p \times \hat{I}_\uparrow \, \hat{I}_\downarrow\, |\rangle =
-1^p \times \hat{I}_\uparrow \, |\rangle \otimes\ \hat{I}_\downarrow \, |\rangle
\]
\[
\mPsi = \sum_\mi \mci |\mi\rangle = \sum_{\mk=1}^{\Ndetup}
\sum_{\mm=1}^{\Ndetdn} C_{\textcolor{darkgreen}{km}} \mDup \mDdn
\]
- If $\mDup$ and $\mDdn$ are represented as $\nmo$ -bit strings, this
transformation can be done in $\mathcal{O}(\Ndet \times \nmo)$.
* Stochastic evaluation of the PT2 correction
** Epstein-Nesbet Second order correction
Consider a wave function $\mPsi$ expanded on an /arbitrary/ set $\mD$
of $\Ndet$ orthonormal Slater determinants,
\[
\mPsi = \sum_{\mi \in \mD} {\mci} | \mi \rangle, \;\;\;
\Evar = \frac{ \langle \mPsi | \mH | \mPsi \rangle}{\langle \mPsi | \mPsi \rangle}
\]
The Epstein-Nesbet 2nd order correction to the energy is
\[
\Ept = \sum_{\ma \in \mA} \frac{ \langle \mPsi | \mH | \ma \rangle
\langle \ma | \mH | \mPsi \rangle}{\Evar - \langle \ma | \mH | \ma \rangle }
\]
The set $\mA$ contains the Slater determinants
- that are not in $\mD$
- for which $d(\mi,\ma) = 1$ or $2$ for at least one pair $(\mi,\ma)$
** Formal scaling
\[
\Ept = \sum_{\ma \in \mA} \frac{ \langle \mPsi | \mH | \ma \rangle
\langle \ma | \mH | \mPsi \rangle}{\Evar - \langle \ma | \mH | \ma
\rangle} = \sum_{\ma \in \mA} \frac{ \left( \sum_{\mi \in \mD} \mci
\langle \mi | \mH | \ma \rangle \right)^2 }{\Evar - \langle \ma |
\mH | \ma \rangle }
\]
- Size of $\mA$ : size of $(\hat{T}_1 + \hat{T}_2)|\mPsi\rangle$
- Number of non-zero terms : $d(\mi,\ma) \le 2 \sim \Ndet \times
\left[ \textcolor{darkgreen}{\left( \nea \times (\nmo - \nea) \right)^2 }\right]$
- *Expensive*
** Solutions to make simulations possible
*** "Non-general'' but /conventional/ solutions:
- Partition the MO space into different classes (active, virtual, inactive, \emph{etc})
- Use another zeroth-order Hamiltonian (CAS-PT2, NEV-PT2)
*** Solutions applicable to /any/ wave function:
- Truncation of $\mD$ to consider only contributions due to large
${\mci}$ \\
But: Truncation $\longrightarrow$ bias because $\Ept$ is a sum of
same-sign values (negative).
- Algorithmic improvement
- Monte Carlo sampling in $\mA$. \\
But: Statistical error decreases as
$\mathcal{O}\left(1/\sqrt{\Nsamples}\right) \Longrightarrow$
Difficult to get $\textcolor{darkgreen}{10^{-5} a.u}$ precision.
- Parallelism
*** :B_ignoreheading:
:PROPERTIES:
:BEAMER_env: ignoreheading
:END:
** Central idea
#+LATEX: \begin{columns}
#+LATEX: \begin{column}{0.5\textwidth}
- Choose an arbitrary ordering of $|\mi \rangle$. Natural choice:
\[w_\mi = \frac{\mci ^2}{\langle \mPsi | \mPsi \rangle } \]
- Make \emph{disjoint} groups $\mA_{\mi}$ of $| \ma \rangle$ originating from the same generator $| \mi \rangle$
- Each $\mA_\mi$ has its own contribution $\epsi$ to $\Ept$
#+LATEX: \end{column}
#+LATEX: \begin{column}{0.5\textwidth}
#+ATTR_LATEX: :width 0.9\columnwidth
[[./fig1.pdf]]
#+LATEX: \end{column}
#+LATEX: \end{columns}
** Central idea
\begin{eqnarray*}
\Ept &=& \sum_{\ma \in \mA} \frac{ \left( \langle \mPsi | \mH | \ma \rangle
\right)^2 }{\Evar - \langle \ma | \mH | \ma \rangle } \\ &=& \sum_{\mi \in \mD}
\sum_{\ma_\mi \in \mA_\mi} \frac{ \left( \langle \mPsi | \mH | \ma_\mi \rangle
\right)^2}{\Evar - \langle \ma_\mi | \mH | \ma_\mi \rangle} \\
&=& \sum_{\mi \in \mD} \epsi
\end{eqnarray*}
*** Contribution per /internal/ determinant
\[
\epsi = \sum_{\ma_\mi \in \mA_\mi} \frac{ \left( \langle \mPsi | \mH | \ma_\mi
\rangle \right)^2}{\Evar - \langle \ma_\mi | \mH | \ma_\mi \rangle}
\]
** From $\mathcal{O}(N_\text{det}^2)$ to $\mathcal{O}(N_\text{det})$
\[
\mPsi = \sum_\mi \mci |\mi\rangle = \sum_{\mk=1}^{\Ndetup} \sum_{\mm=1}^{\Ndetdn} C_{\textcolor{red}{km}} \mDup \mDdn
\]
- Sorting is $\mathcal{O}(\Ndet)$
- $\langle \mi | \mH | \ma \rangle \langle \ma | \mH | \mj \rangle = 0$ when $d(\mi, \mj) > 4$
- Loop over $\Ndetup$ determinants (rows of the $C$ matrix) \\
Remove all the rows where $d(\mDup, \mDupp)>4$ ($\sim \mathcal{O}(\sqrt{\Ndet})$)
- Loop over $\Ndetdn$ determinants (columns of the $C$ matrix) \\
Remove all the columns where $d(\mDdn, \mDdnn)>4)$
- The remaining number of determinants is bounded by the size of the CISDTQ space
** Finer filtering
\[
\epsi = \sum_{\ma \in \mA_\mi} \frac{ \langle \mPsi | \mH | \ma_\mi \rangle
\langle \ma_\mi | \mH | \mPsi' \rangle}{\Evar - \langle \ma_\mi | \mH | \ma_\mi
\rangle}
\]
- We know that all the $| \ma_\mi \rangle$ are singles and doubles
with respect to $| \mi \rangle$
- $|\mPsi'\rangle$ is the projection of $\mPsi$ on the subspace of
determinants in $\mD$ which are no more than quadruply excited
with respect to $| \mi \rangle$
- For a subset of excitations $\mi \rightarrow \ma$,
$|\mPsi'\rangle$ is filtered further with possible hole/particle constraints
** Monte Carlo sampling
\[
\epsi = \sum_{\ma_\mi \in \mA_\mi} \frac{ \left( \langle \mPsi | \mH | \ma_\mi
\rangle \right)^2}{\Evar - \langle \ma_\mi | \mH | \ma_\mi \rangle}
\]
1. $\langle \mPsi | \mH | \ma_\mi \rangle = \sum_{\mj \ge \mi} \mcj\, \langle\, \mj | \mH | \ma_\mi \rangle$
2. $\langle \ma_\mi | \mH | \ma_\mi \rangle$ is always large
(otherwise $|\ma_\mi \rangle$ would be better in the
\textcolor{red}{variational space}, and PT is questionable)
*** :B_ignoreheading:
:PROPERTIES:
:BEAMER_env: ignoreheading
:END:
- $\forall \mi \in \mD : \epsi \le 0$
- $|\epsi|$ is expected to decrease as $\mci^2$
- The computational cost decreases with $\mi$
*** Monte Carlo formulation
\[
\Ept = \sum_{\mi \in \mD} \epsi = \sum_{\mi \in \mD} \mpi \frac{\epsi}{\mpi} =
\left\langle \frac{\epsi}{\mpi} \right\rangle_{\mpi}
\]
** Naive sampling
#+LATEX: \begin{columns}
#+LATEX: \begin{column}{0.2\textwidth}
Uniform sampling:
\mpi = \frac{1}{\Ndet}$
#+LATEX: \end{column}
#+LATEX: \begin{column}{0.8\textwidth}
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./dist2_noise.png]]
#+LATEX: \end{column}
#+LATEX: \end{columns}
** Improved sampling
#+LATEX: \begin{columns}
#+LATEX: \begin{column}{0.2\textwidth}
Sampling : $\mpi = \mci^2$
#+LATEX: \end{column}
#+LATEX: \begin{column}{0.8\textwidth}
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./eici2.png]]
#+LATEX: \end{column}
#+LATEX: \end{columns}
** Lazy evaluation
Only $\Ndet$ contributions $\epsi \longrightarrow$ all $\epsi$ can be stored in memory.
*** Lazy Evaluation (Wikipedia)
In programming language theory, /lazy evaluation/, or /call-by-need/ is
an evaluation strategy which delays the evaluation of an expression until its
value is needed (non-strict evaluation) and which also avoids repeated
evaluations (sharing).
*** :B_ignoreheading:
:PROPERTIES:
:BEAMER_env: ignoreheading
:END:
#+BEGIN_SRC python
def lazy_e(i):
if not e_is_computed[i]:
e[i] = compute_e(i)
e_is_computed[i] = true
return e[i]
#+END_SRC
** Monte Carlo with Lazy Evaluation
\[
\Ept = \sum_{\mi \in \mD} \epsi = \sum_{\mi \in \mD} \mpi \frac{\epsi}{\mpi} =
\left\langle \frac{\epsi}{\mpi} \right\rangle_{\mpi}
\]
- Draw a generator determinant $|\mi\,\rangle$ with probability $\mpi$
- Increment $n_\mi$, the number of evaluations of $\epsi$
- If $\epsi$ is not already computed, compute it and store its value
- $\Ept \sim \sum_{\mi \in \mD} \frac{n_\mi}{\Nsamples} \frac{\epsi}{\mpi}$
- Statistical error : $\mathcal{O}\left(1/\sqrt{\Nsamples}\right)$
- Lazy evaluation : Exponential acceleration (time to solution)
** Monte Carlo with Lazy Evaluation
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./samples.pdf]]
** Monte Carlo with Lazy Evaluation
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./lazy_e.pdf]]
** Monte Carlo with Lazy Evaluation
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./lazy_err.pdf]]
** Monte Carlo with Variance reduction
- Noise can be smoothed out by averaging
- Split $\mD$ into $\mM$ \emph{equiprobable} sets : "Comb"
\[
\Ept = \sum_{\mi \in \mD} \epsi = \sum_{\mk=1}^{\mM} \sum_{\mi_\mk \in \mD_\mk} \epsik
\]
*** New Monte Carlo estimator
\[
\Ept = \left \langle
\frac{1}{\mM} \sum_{\mk=1}^{\mM} \frac{\epsik}{{\textcolor{red}{p_{I_k}}}}
\right \rangle_{\textcolor{red}{(p_{I_1}, \dots, p_{I_M})}}
\]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./Comb1.pdf]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./Comb2.pdf]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./Comb3.pdf]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./dist2.png]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./lazy_e.pdf]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./comb_e.pdf]]
** Monte Carlo with Variance reduction
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./comb_err.pdf]]
** Hybrid deterministic/stochastic scheme
- When all the determinants have been drawn, the \emph{exact} $\Ept$ can be computed
- $\Longrightarrow$ The result with zero statistical error can be reached
in a finite time
- In typical wave functions, $90\%$ of the norm is on a few determinants
- Compute the few first contributions $\epsi$, and perform the MC in the rest
\[
\Ept = \sum_{\mi \in \mD_D} \epsi + \left \langle
\frac{1}{\mM} \sum_{\mk=1}^{\mM} \frac{\epsik}{\textcolor{red}{p_{I_k}}}
\right \rangle_{(\textcolor{red}{p_{I \in \mD_S})}}
\]
** Hybrid deterministic/stochastic scheme
Make the deterministic part grow during the calculation.
*** At each MC step
- Draw a random number
- Find the determinants selected by the comb (increment $n_\mi$'s)
- Compute the $\epsi$ which have not been yet computed
- Compute deterministically the first non-computed determinant
- If a tooth of the comb is completely filled $\Longrightarrow$ Deterministic
*** At any time
\[
\Ept(t) = \sum_{\mi \in \mD_D(t)} \epsi + \sum_{\mi \in \mD_S(t)} \frac{1}{\mM(t)}
\frac{n_\mi(t)}{\Nsamples(t)} \frac{\epsi}{\mpi}
\]
** Hybrid deterministic/stochastic scheme
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./lazy_e.pdf]]
** Hybrid deterministic/stochastic scheme
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./comb_e.pdf]]
** Hybrid deterministic/stochastic scheme
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./hybrid_e.pdf]]
** Hybrid deterministic/stochastic scheme
#+CAPTION: F$_2$, cc-pVDZ, \textcolor{red}{$10^6$} determinants in the variational space
#+ATTR_LATEX: :height 0.75\textheight
[[./hybrid_err.pdf]]
** Some timings: Cr$_2$, $2\,10^7$ determinants, 800 cores
|---------+-------------------------------------+--------------------------|
| Basis | $\Ept$ | Wall-clock time |
|---------+-------------------------------------+--------------------------|
| cc-pVDZ | \textcolor{red}{$-0.068\,3(1)$} | 14 min |
| | \textcolor{red}{$-0.068\,36(1)$} | 55 min |
| | \textcolor{red}{$-0.068\,361(1)$} | 2.4 hr |
| | \textcolor{red}{$-0.068\,360\,604$} | 3 hr |
|---------+-------------------------------------+--------------------------|
| cc-pVTZ | \textcolor{red}{$-0.124\,4(5)$} | 19 min |
| | \textcolor{red}{$-0.124\,7(1)$} | 58 min |
| | \textcolor{red}{$-0.124\,63(1)$} | 3.5 hr |
| | \textcolor{red}{$-0.124\,642(1)$} | 8.7 hr |
| | --- | $\sim$ 15 hr (estimated) |
|---------+-------------------------------------+--------------------------|
| cc-pVQZ | \textcolor{red}{$-0.155\,8(5)$} | 56 min |
| | \textcolor{red}{$-0.155\,9(1)$} | 2.5 hr |
| | \textcolor{red}{$-0.155\,95(1)$} | 9.0 hr |
| | \textcolor{red}{$-0.155\,952(1)$} | 18.5 hr |
| | --- | $\sim$ 29 hr (estimated) |
|---------+-------------------------------------+--------------------------|
** Parallel efficiency
#+ATTR_LATEX: :height 0.8\textheight
[[./speedup_pt2.png]]
* Export :noexport:
#+BEGIN_SRC elisp :output none
(setq org-latex-listings 'minted
org-latex-packages-alist '(("" "minted"))
org-latex-pdf-process
'("pdflatex -shell-escape -interaction nonstopmode -output-directory %o %f"
"pdflatex -shell-escape -interaction nonstopmode -output-directory %o %f"
"pdflatex -shell-escape -interaction nonstopmode -output-directory %o %f"))
(setq org-latex-minted-options '(("breaklines" "true")
("breakanywhere" "true")))
(setq org-latex-minted-options
'(("frame" "lines")
("fontsize" "\\scriptsize")
("linenos" "")))
(org-beamer-export-to-pdf)
#+END_SRC
#+RESULTS:
: /home/scemama/TEX/Nuclearistes/algorithms.pdf