\section{Transport and depletion}
\subsection{Particle transport}
The particle transport describes the evolution of the particles concentration in a certain position in the phase space $(\Vec{r},E,\Vec{\Omega})$, by considering all the mechanisms of production, depletion and transport inside the domain \cite{Bell1970NuclearReactorTheory}. The governing equation for neutron flux is presented in eq. \eqref{transport}.
\begin{multline}
   \frac{1}{v} \frac{\partial \phi (\Vec{r},E,\Vec{\Omega},t)}{\partial t} \ + \ \Vec{\Omega} \; \nabla \phi (\Vec{r},E,\Vec{\Omega},t) \ + \  \sigma_t (E) N_j (\Vec{r},t) \;\phi (\Vec{r},E,\Vec{\Omega},t) \ = \\ 
    \int d\Vec{\Omega'} \int dE'\; \sigma_s (E) N_j (\Vec{r},t) \;f_s (\Vec{r},E' \rightarrow E, \Vec{\Omega} \cdot \Vec{\Omega'})\; \phi (\Vec{r},E,\Vec{\Omega},t) \ + \  Q (\Vec{r},E,\Vec{\Omega},t);
    \label{transport}
\end{multline}

The variable that appears in the equation is the neutron angular flux, defined as: $$\phi (\Vec{r},E,\Vec{\Omega},t) \ = \ v \; N (\Vec{r},E,\Vec{\Omega},t). $$ 
Three independent variables define the phase space: $\Vec{r}$ is the position, E is the energy and $\Vec{\Omega}$ is the solid angle $\Vec{v} = v \, \Vec{\Omega}$ determining the direction of the motion.
The meaning of the microscopic and macroscopic has already been discussed and $N_j$ is the concentration of the background.
$f_s$ is the scattering probability density function and it represents the probability that a neutron of direction $\Vec{\Omega'}$ and energy $E'$ has to scatter in the direction interval $d\Vec{\Omega}$ about $\Vec{\Omega}$ with energy in $dE$ about $E$.
$Q$ represents the source term, that can be of different nature; fixed source providing particles from outside the system and fission processes producing neutrons if fissile nuclides are present in the system.
The number of neutrons produced by fission is not constant and it is a function of the target nuclide and of the energy, $\nu (\Vec{r},E')$.
The energy of produced neutrons is determined by the probability function $\chi (\Vec{r},E)$ and the process is isotropic so the probability density function for the angular  distribution is 
\begin{equation}
  \int d\Vec{\Omega} \ f(\Vec{\Omega}) = 1 \ \rightarrow \ f(\Vec{\Omega}) \ = \ \frac{1}{4 \pi}.  
    \label{angular_distr}
\end{equation}\\

Therefore, the fission term appearing in the transport equation can be written as
\begin{equation}
   \int d\Vec{\Omega'} \int dE' \ \frac{\chi (\Vec{r},E)}{4 \pi} \ \nu (\Vec{r},E') \ \sigma_f(E) N_i(\Vec{r},t) \ \phi (\Vec{r},E,\Vec{\Omega},t).
\end{equation}

In a critical system the leakages are fully compensated by neutrons produced by fission and the reduced mathematical formulation is the following:
\begin{equation}
    \label{tr_crit}
    \hat{L} \, \phi (\Vec{r},E,\Vec{\Omega}) = \frac{1}{k} \, \hat{F} \, \phi (\Vec{r},E,\Vec{\Omega}).
\end{equation}
In the formulation of eq. \eqref{tr_crit} the term $k$, called multiplication factor \cite{Lamarsh} is introduced; it is the eigenvalue of the equation and associated eigenvector is the solution of the problem $\phi$. This mathematical trick allows to make the system critical just adjusting the number of neutrons produced by fission.\\
Sub-critical system have $k\,<\,1$, the geometrical configuration does not allow criticality and an external source fills the gap.
This is the case of an \acrshort{ADS}, like MYRRHA, in which the surplus of neutrons if provided by the spallation source.

\subsection{Depletion equations}
Bateman equations \cite{Bateman2011SolutionTransformations} describe the evolution of nuclide concentrations in a system undergoing irradiation and radioactive decay. They are a generalization of the exponential decay law and they include every reaction to assess the gross production and depletion of a certain nuclide.
The general mathematical formulation, taking into account different particle fluxes and spatial dependence is presented in eq. \eqref{bateman_full}.
\begin{multline}
    \frac{\partial N_i (\Vec{r},t)}{\partial t} \ = \ \sum_j (\lambda_{j \rightarrow i} + \sum_M \int dE \sigma_{j \rightarrow i}^M (E) \phi^M(\Vec{r},E)) N_j (\Vec{r},t) - \\
    (\lambda_{i} + \sum_M \int dE \sigma_{tot,i}^M (E) \phi^M(\Vec{r},E)) N_i (\Vec{r},t)
\label{bateman_full}
\end{multline}
The following notation is used:
\begin{itemize}
	\item $i$, $j$ are the indices of nuclides;
	\item $\lambda_{j \rightarrow i}$  is the decay constant in which $j$ is the mother and $i$ the daughter;
	\item $\lambda_{i}$  is the decay constant of $i$;
	\item $M$ stands for the particle fluxes involved;
    \item $\sigma_{j \rightarrow i}^M$ is the microscopic cross section for the reaction producing $i$ when $j$ interact with the particle $M$;
    \item $\sigma_{tot,i}^M$ is the removal cross section of $i$ when interacting with $M$;
    \item $\phi^M(\Vec{r},E))$ is the flux of particle $M$.
\end{itemize}
The closed set of first order integro-differential equations can be written in matrix form and solved by numerical scheme, its solution is the evolution of the nuclide vector in the system ($\Vec{N} (t)$).
\subsection{Coupling transport and depletion}
The variables $N_j(\Vec{r},t)$ in \ref{transport} and $\phi (\Vec{r},E,t)$ in \eqref{bateman_full} establish the coupling between transport and depletion since the spectrum coming from transport simulation is necessary to estimate the reaction rates appearing in the Bateman equations

\begin{equation} \label{rr}
 <\sigma_{j \rightarrow i} \phi > = \int dE \sigma_{j \rightarrow i} (E) \phi(\Vec{r},E)
\end{equation}
and the updated nuclide vector must be fed to the transport model and it affects the new flux.\\
The time coupling is therefore often solved by separation of spatial and temporal components of the solution. This approach evaluates both the neutronic evolution and the nuclides concentrations in the system.
\subsection{Derived quantities}
The nuclide vector is employed to estimate other quantities of major importance for the analysis such as activity and decay heat. The activity is defined as the number of decays per unit time so it is proportional to the atomic density, multiplied for the decay constant $\lambda$. So the activity evolution of a certain nuclide mimics the concentration evolution. The same applies for the decay heat rate, the other quantity relevant in this work. It consists of the energy released in the decay $Q_i$ multiplied by the decay rate of a nuclide (that is the activity).\\
Mathematical expressions of the two quantities are shown in Eqs. \eqref{act_dh}.
\begin{equation}
    \begin{gathered}
    A(t) = \sum_i \lambda_i\,N_i(t)\\
    DH(t) = \sum_i Q_i\,\lambda_i\,N_i(t)
    \label{act_dh}
    \end{gathered}
\end{equation}
