source: trunk/documents/UserDoc/UsersGuides/PhysicsReferenceManual/latex/hadronic/theory_driven/Evaporation/ModelDescription.tex @ 1208

\section{Model description.}

\hspace{1.0em} The Weisskopf treatment is an application of the
detailed balance principle that relates the probabilities to go from
a state $i$ to another $d$ and viceversa through the density of states
in the two systems:
%%
\begin{equation}
  \label{evap:1}
  P_{i \rightarrow d} \rho(i) = P_{d \rightarrow i} \rho(d)
\end{equation}
%%
where $P_{d \rightarrow i}$ is the probability per unit of time of a
nucleus $d$ captures a particle $j$ and form a compound nucleus $i$
which is proportional to the compound nucleus cross section
$\sigma_{\mathrm{inv}}$. Thus, the probability that a parent nucleus
$i$ with an excitation energy $E^*$ emits a particle $j$ in its ground
state with kinetic energy $\varepsilon$ is 
%%
\begin{equation}
  \label{evap:2}
  P_j(\varepsilon) \mathrm{d}\varepsilon = g_j
  \sigma_{\mathrm{inv}}(\varepsilon) \frac{\rho_d(E_{\mathrm{max}} -
  \varepsilon)}{\rho_i(E^*)} \varepsilon \mathrm{d}\varepsilon
\end{equation}
%%
where $\rho_i(E^*)$ is the level density of the evaporating nucleus,
$\rho_d(E_{\mathrm{max}} - \varepsilon)$ that of the daugther
(residual) nucleus after emission of a fragment $j$ and
$E_{\mathrm{max}}$ is the maximum energy that can be carried by the
ejectile. With the spin $s_j$ and the mass $m_j$ of the emitted
particle, $g_j$ is expressed as $g_j = ( 2 s_j + 1 ) m_j / \pi^2
\hbar^2$.

This formula must be implemented with a suitable form for the level
density and inverse reaction cross section. We have followed, like
many other implementations, the original work of Dostrovsky \textit{et al.}
 \cite{evap.Dostrovsky59} (which represents the first Monte Carlo
code for the evaporation process) with slight modifications. The
advantage of the Dostrovsky model is that it leds to a simple
expression for equation \ref{evap:2} that can be analytically
integrated and used for Monte Carlo sampling.


\subsection{Cross sections for inverse reactions.}

\hspace{1.0em} The cross section for inverse reaction is expressed by
means of empirical equation \cite{evap.Dostrovsky59}
%%
\begin{equation}
  \label{evap:3}
  \sigma_{\mathrm{inv}}(\varepsilon) = \sigma_g \alpha \left( 1 + 
  \frac{\beta}{\varepsilon} \right)  
\end{equation}
%%
where $\sigma_g = \pi R^2$ is the geometric cross section. 

In the case of neutrons, $\alpha = 0.76+2.2A^{-\frac{1}{3}}$ and
$\beta = (2.12 A^{-\frac{2}{3}} - 0.050)/\alpha$ MeV. This equation
gives a good agreement to those calculated from continuum theory
\cite{evap.Blatt52} for intermediate nuclei down to $\varepsilon \sim 0.05$
MeV. For lower energies $\sigma_{\mathrm{inv},n}(\varepsilon)$ tends
toward infinity, but this causes no difficulty because only the
product $\sigma_{\mathrm{inv},n}(\varepsilon)\varepsilon$ enters in
equation \ref{evap:2}. It should be noted, that the inverse cross
section needed in  \ref{evap:2} is that between a neutron of kinetic
energy $\varepsilon$ and a nucleus in an excited state. 

For charged particles (p, d, t, $^3$He and $\alpha$), $\alpha =
(1+c_j)$ and $\beta = -V_j$, where $c_j$ is a set of parameters
calculated by Shapiro \cite{evap.Shapiro53} in order to provide a good fit
to the continuum theory \cite{evap.Blatt52} cross sections and $V_j$ is the
Coulomb barrier.

\subsection{Coulomb barriers.}

\hspace{1.0em} Coulomb repulsion, as calculated from elementary
electrostatics are not directly applicable to the computation of
reaction barriers but must be corrected in several ways. The first
correction is for the quantum mechanical phenomenoon of barrier
penetration. The proper quantum mechanical expressions for barrier
penetration are far too complex to be used if one wishes to retain
equation \ref{evap:2} in an integrable form. This can be approximately
taken into account by multiplying the electrostatic Coulomb barrier by
a coefficient $k_j$ designed to reproduce the barrier penetration
approximately whose values are tabulated \cite{evap.Shapiro53}.
%%
\begin{equation}
  \label{evap:4}
  V_j = k_j \frac{Z_j Z_d e^2}{R_c}
\end{equation}
%%
The second correction is for the separation of the centers of the
nuclei at contact, $R_c$. We have computed this separation as $R_c =
R_j + R_d$ where $R_{j,d} = r_c A_{j,d}^{1/3}$ and $r_c$ is
given \cite{evap.Iljinov94} by
%%
\begin{equation}
  \label{evap:5}
  r_c = 2.173 \frac{1 + 0.006103 Z_j Z_d}{1 + 0.009443 Z_j Z_d}
\end{equation}
%%
 

\subsection{Level densities.}

\hspace{1.0em} The simplest and most widely used level density based
on the Fermi gas model  are those of Weisskopf \cite{evap.Weisskopf37} for
a completely degenerate Fermi gas. We use this approach with the
corrections for nucleon pairing proposed by Hurwitz and Bethe
\cite{evap.Hurwitz51} which takes into account the displacements of the
ground state:
%%
\begin{equation}
  \label{evap:6}
  \rho(E) = C \exp{\left( 2 \sqrt{a(E-\delta)} \right)}
\end{equation}
%%
where $C$ is considered as constant and does not need to be specified
since only ratios of level densities enter in equation \ref{evap:2}.
$\delta$ is the pairing energy correction of the daughter nucleus
evaluated by Cook \textit{et al.} \cite{evap.Cook67} and Gilbert and
Cameron \cite{evap.Gilbert65} for those values not evaluated by Cook
\textit{et al.}. The level density parameter is calculated according
to:
%%
\begin{equation}
  \label{evap:7}
  a(E,A,Z) = \tilde{a}(A) \left \{ 1 + \frac{\delta}{E} [1 - \exp(-\gamma
  E)] \right\}
\end{equation}
%%
and the parameters calculated by Iljinov \textit{et al.}
\cite{evap.Iljinov92} and shell corrections of Truran, Cameron and Hilf
\cite{evap.Truran70}.

\subsection{Maximum energy available for evaporation.}

\hspace{1.0em} The maximum energy avilable for the evaporation process
(\textit{i.e.} the maximum kinetic energy of the outgoing fragment) is
usually computed like $E^* - \delta - Q_j$ where is the separation
energy of the fragment $j$: $Q_j = M_i - M_d - M_j$ and $M_i$, $M_d$
and $M_j$ are the nclear masses of the compound, residual and
evporated nuclei respectively. However, that expression does not
consider the recoil energy of the residual nucleus. In order to take
into account the recoil energy we use the expression
%%
\begin{equation}
  \label{evap:8}
  \varepsilon_j^{\mathrm{max}} = \frac{(M_i + E^* - \delta)^2 + M_j^2
  - M_d^2}{2 (M_i + E^* - \delta)} - M_j
\end{equation}

\subsection{Total decay width.}

\hspace{1.0em} The total decay width for evaporation of a fragment $j$
can be obtained by integrating equation \ref{evap:2} over kinetic
energy
%%
\begin{equation}
  \label{evap:9}
  \Gamma_j = \hbar \int_{V_j}^{\varepsilon_j^{\mathrm{max}}}
  P(\varepsilon_j) \mathrm{d}\varepsilon_j
\end{equation}
%%
This integration can be performed analiticaly if we use equation
\ref{evap:6} for level densities and equation \ref{evap:3} for inverse
reaction cross section. Thus, the total width is given by
%%
\begin{eqnarray}
\label{eq:10}
  \Gamma_j  = \frac{g_j m_j R_d^2}{2 \pi \hbar^2} \frac{\alpha}{a_d^2}
    & \times &
    \Biggl \lgroup \biggl \{ \left(\beta a_d - \frac{3}{2}\right) + a_d
    (\varepsilon_j^{\mathrm{max}} - V_j) \biggr \}
    \exp{\left\{-\sqrt{a_i(E^*-\delta_i)}\right\}} + \nonumber \\
    & & \biggl \{ (2\beta a_d - 3) \sqrt{a_d
       (\varepsilon_j^{\mathrm{max}} - V_j)} + 2 a_d
    (\varepsilon_j^{\mathrm{max}} - V_j) \biggr \} \times \nonumber \\
    & & \exp{
       \left\{
       2 \left[
         \sqrt{a_d(\varepsilon_j^{\mathrm{max}} - V_j)}
         -  
         \sqrt{a_i(E^* - \delta_i)}
       \right]
       \right\}
       }
     \Biggr \rgroup
\end{eqnarray}
%%
where $a_d = a(A_d,Z_d,\varepsilon_j^{\mathrm{max}})$ and $a_i =
a(A_i,Z_i,E^*)$.

\section{GEM Model}

\hspace{1.0em} As an alternative model we have implemented the
generalized evaporation model (GEM) by Furihata
\cite{evap.Furihata00}. This model considers emission of fragments heavier
than $\alpha$ particles and uses a more accurate level density
function for total decay width instead of the approximation used by
Dostrovsky. We use the same set of parameters but for heavy ejectiles
the parameters determined by Matsuse \textit{et al.} \cite{evap.Matsuse82}
are used.

Based on the Fermi gas model, the level density function is expressed
as
%%
\begin{equation}
\label{evap:11}
  \rho(E) =
    \left\{ \begin{array}{ll}
      \frac{\sqrt{\pi}}{12} 
      \frac{e^{2\sqrt{a(E-\delta)}}}{a^{1/4}(E-\delta)^{5/4}} &  
      \rm{for} \quad E \geq E_x \\
      \frac{1}{T} e^{(E-E_0)/T} & \rm{for} \quad E < E_x
    \end{array} \right.
\end{equation}
%%
where $E_x = U_x + \delta$ and $U_x = 150/M_d + 2.5$ ($M_d$ is the
mass of the daughter nucleus). Nuclear temperature $T$ is given as
$1/T = \sqrt{a/U_x} - 1.5U_x$, and $E_0$ is defined as $E_0 = E_x -
T(\log{T} - \log{a}/4 - (5/4)\log{U_x} + 2 \sqrt{aU_x})$.

By substituting equation \ref{evap:11} into equation \ref{evap:2} and
integrating over kinetic energy can be obtained the following
expression
%%
\begin{equation}
\label{evap:12}
  \Gamma_j = \frac{ \sqrt{\pi} g_j \pi R_d^2 \alpha}{12 \rho(E^*)} \times
  \left\{ \begin{array}{ll}
    \{I_1(t,t) + (\beta+V)I_0(t)\} &
    \rm{for} \quad \varepsilon_j^{\mathrm{max}} - V_j < E_x \\
    \{I_1(t,t_x)+I_3(s,s_x)e^s+ & \\
    (\beta+V)(I_0(t_x)+I_2(s,s_x)e^s)\} &
    \rm{for} \quad \varepsilon_j^{\mathrm{max}} - V_j \geq E_x .
  \end{array} \right.
\end{equation}
%%
$I_0(t)$, $I_1(t,t_x)$, $I_2(s,s_x)$, and $I_3(s,s_x)$ are expressed
as:
\begin{eqnarray}
  I_0(t)     & = & e^{-E_0/T} (e^t -1) \\
  I_1(t,t_x) & = & e^{-E_0/T} T \{(t - t_x + 1)e^{t_x} - t -1 \} \\ 
  I_2(s,s_x) & = & 2\sqrt{2} \biggl \{ s^{-3/2} + 1.5  s^{-5/2} + 3.75
  s^{-7/2} - \nonumber \\ 
             &   & (s_x^{-3/2} + 1.5 s_x^{-5/2} + 3.75 s_x^{-7/2}) \biggr \} \\
  I_3(s,s_x) & = & \frac{1}{2\sqrt{2}} \Biggl [ 2 s^{-1/2} + 4
                  s^{-3/2} + 13.5  s^{-5/2} + 60.0  s^{-7/2} + \nonumber \\ 
             &   & 325.125 s^{-9/2} - 
  \biggl \{ 
  (s^2 - s_x^2) s_x^{-3/2} + 
  (1.5s^2 + 0.5s_x^2) s_x^{-5/2} + \nonumber \\
  & & 
  (3.75s^2 + 0.25s_x^2) s_x^{-7/2} + 
  (12.875s^2 + 0.625s_x^2) s_x^{-9/2} + \nonumber \\ 
  & & 
  (59.0625s^2 + 0.9375s_x^2) s_x^{-11/2} + \nonumber \\
  & &
  (324.8s^2 + 3.28s_x^2) s_x^{-13/2} +
  \biggr \}
\Biggr ]
\end{eqnarray}
%%
where $t = (\varepsilon_j^{\mathrm{max}}-V_j)/T$, $t_x = E_x/T$, $s = 2
\sqrt{a(\varepsilon_j^{\mathrm{max}}-V_j -\delta_j)}$ and $s_x = 2
\sqrt{a(E_x - \delta)}$.

Besides light fragments, 60 nuclides up to $^{28}$Mg are considered,
not only in their ground states but also in their exited states, are
considered. The excited state is assumed to survive if its lifetime
$T_{1/2}$ is longer than the decay time, \textit{i. e.},
$T_{1/2}/\ln{2} > \hbar/\Gamma_j^*$, where $\Gamma_j^*$ is the
emission width of the resonance calculated in the same manner as for
ground state particle emission. The total emission width of an
ejectile $j$ is summed over its ground state and all its excited
states which satisfy the above condition.

%%% Local Variables: 
%%% mode: latex
%%% TeX-master: "EvaporationModel"
%%% End: