diff --git a/.gitignore b/.gitignore
index 90400f6c76aa63f687ba26ad224bb5a0735398a6..2f77a0a9ab29ef7ca65377fb518273ada7389d7e 100644
--- a/.gitignore
+++ b/.gitignore
@@ -16,4 +16,5 @@ sammy/TriBITS
.DS_Store
sammy/[bB][uU][iI][lL][dD]*
/install
+[Ss][Cc][Aa][Ll][Ee]
.vscode
diff --git a/docs/tex/sammy_main.tex b/docs/tex/sammy_main.tex
index da8f82f4596f04e44d2036efdf02a31a0308f643..8a7b1c8312fb4f79a2f51137dee3b0663bf9312e 100644
--- a/docs/tex/sammy_main.tex
+++ b/docs/tex/sammy_main.tex
@@ -50,6 +50,7 @@
\usepackage{lipsum}
\usepackage[hidelinks]{hyperref} % ,backref
\usepackage[all]{nowidow}
+\usepackage{threeparttable}
%\DeclareMathAlphabet{\pazocal}{OMS}{zplm}{m}{n}
\newcommand{\Ma}{\mathcal{M}}
\newcommand{\Mb}{\pazocal{M}}
diff --git a/docs/tex/sammyrefs.bib b/docs/tex/sammyrefs.bib
index aca4cfa1ebcb5b25275756870710358fa7722189..fdd0f044f96f4c1c3caa2336b7c82e25ad9a4bbd 100644
--- a/docs/tex/sammyrefs.bib
+++ b/docs/tex/sammyrefs.bib
@@ -1161,7 +1161,6 @@ abstract = "The response of the RPI 16-segment NaI(Tl) multiplicity detector sys
number = {CONF-880546},
}
-<<<<<<< HEAD
%----------- SAMMY specific references ---------------
@conference{larson_multscat_2002,
@@ -1191,9 +1190,12 @@ abstract = "The response of the RPI 16-segment NaI(Tl) multiplicity detector sys
publisher = {Los Alamos Scientific Laboratory},
location = {Los Alamos, New Mexico, USA},
}
-=======
-
->>>>>>> f90b3e53... Added tex files for the manual multiple scattering documentation
+@article{lynn_1968,
+ title={The theory of neutron resonance reactions},
+ author={Lynn, John Eric},
+ year={1968},
+ publisher={Clarendon Press, Oxford}
+}
@Comment{jabref-meta: databaseType:bibtex;}
diff --git a/docs/tex/scattering-theory.tex b/docs/tex/scattering-theory.tex
index d1d27cc67a8eb55ef4959c9fbb9575c0fa58a459..8bc8d065e1a71d88f7d5090c8f8fc09c2ccf6116 100644
--- a/docs/tex/scattering-theory.tex
+++ b/docs/tex/scattering-theory.tex
@@ -26,26 +26,16 @@ R-matrix theory is expressed in terms of channels, where a channel is defined as
\begin{figure}
\centering
- \includegraphics[width=0.5\textwidth]{figures/rmatrix_inc_channel.pdf} \\[\smallskipamount]
- \includegraphics[width=0.5\textwidth]{figures/rmatrix_interior_region.pdf}\hfill
- \includegraphics[width=0.5\textwidth]{figures/rmatrix_exit_channel.pdf}
+ \includegraphics[width=0.4\textwidth]{figures/rmatrix_inc_channel.pdf} \\[\smallskipamount]
+ \includegraphics[width=0.4\textwidth]{figures/rmatrix_interior_region.pdf}\hfill
+ \includegraphics[width=0.4\textwidth]{figures/rmatrix_exit_channel.pdf}
\caption{Schematic of entrance and exit channels as used in scattering theory. For the interior region (with separation distance $r < a$), no assumptions are made about the nature of the interaction. In the figure, $m$, $i$, and $z$ refer to the mass, spin, and charge of the incident particle while $M$, $I$ and $Z$ refer to the target particle. Orbital angular momentum is denoted by $l$ and velocity by $v$. Primes are used for post-collision quantities.}
\label{scattering-theory_rmatrix_channel_diagram}
\end{figure}
-In Section II.A\ref{}, general equations of scattering theory are presented and their derivations
-discussed. The fundamental R-matrix equations are presented. Section II.A.1\ref{} gives a detailed
-derivation of the equations for a simple case. Section II.A.2\ref{} shows the relationship between the
-R-matrix and the A-matrix, which is another common representation of scattering theory.
+In Section \ref{sec:equations-for-scattering-theory}, general equations of scattering theory are presented and their derivations discussed. The fundamental R-matrix equations are presented. Section II.A.1\ref{} gives a detailed derivation of the equations for a simple case. Section II.A.2\ref{} shows the relationship between the R-matrix and the A-matrix, which is another common representation of scattering theory.
-The approximations to R-matrix theory available in the SAMMY code are detailed in
-Section II.B\ref{}. The recommended choice for most applications is the Reich-Moore approximation,
-described in Section II.B.1\ref{}. For some applications, the Reich-Moore approximation is
-inadequate; for those cases, a method for using SAMMY's Reich-Moore approximation to
-mimic the full (exact) R-matrix is presented Section II.B.2\ref{}. Two historically useful but now
-obsolete approximations are single-level and multilevel Breit Wigner (SLBW and MLBW),
-discussed in Section II.B.3\ref{}. Provisions for including non-compound (direct) effects are
-discussed in Section II.B.4\ref{}.
+The approximations to R-matrix theory available in the SAMMY code are detailed in Section II.B\ref{}. The recommended choice for most applications is the Reich-Moore approximation, described in Section II.B.1\ref{}. For some applications, the Reich-Moore approximation is inadequate; for those cases, a method for using SAMMY's Reich-Moore approximation to mimic the full (exact) R-matrix is presented Section II.B.2\ref{}. Two historically useful but now obsolete approximations are single-level and multilevel Breit Wigner (SLBW and MLBW), discussed in Section II.B.3\ref{}. Provisions for including non-compound (direct) effects are discussed in Section II.B.4\ref{}.
In Section II.C\ref{}, details are given for the SAMMY nomenclature and other conventions,
for transformations to the center-of-momentum system, and for the calculation of penetrability,
@@ -75,6 +65,113 @@ vector sum of the spins of the two particles of the pair: $\vec{s} = \vec{i} + \
of $l$ and $s$: $\vec{J} = \vec{l} + \vec{s}$.
\end{itemize}
+Only $J$ and its associated parity $\pi$ are conserved for any given interaction. The other quantum numbers may differ from channel to channel, as long as the sum rules for spin and parity are obeyed. Within this document and within the SAMMY code, the set of all channels with the same $J$ and $\pi$ are called a ``spin group.''
+
+In all formulae given below, spin quantum numbers (e.g., $J$ ) are implicitly assumed to include the associated parity. Quantized vector sum rules are implicitly assumed to be obeyed. Readers unfamiliar with these sum rules are referred to Section II.C.1.a\ref{} for a mini-tutorial on the subject.
+
+Let the angle-integrated cross sections from entrance channel $c$ to exit channel $c'$ with total angular momentum $J$ be represented by $\sigma_{cc'}$. This cross section is given in terms of the scattering matrix $U_{cc'}$ as
+
+\begin{equation}\label{eq:sigma-ccprime}
+\sigma_{cc'} = \frac{\pi}{k_\alpha^2}g_{J\alpha}\left|e^{2iw_c}\delta_{cc'}-U_{cc'}\right|^2\delta_{JJ'} \:,
+\end{equation}
+
+\noindent
+where $k_\alpha$ is the wave number (and $K_\alpha = \hbar k_\alpha =$ center-of-mass momentum) associated with incident particle pair $\alpha$, $g_{J\alpha}$ is the spin statistical factor, and $w_c$ is the Coulomb phase-shift difference. Note that $w_c$ is zero for non-Coulomb channels. (Details for the charged-particle case are presented in Section II.C.4.\ref{}) The spin statistical factor $g_{J\alpha}$ is given by
+
+\begin{equation}\label{eq:spin-stat-factor}
+g_{J\alpha} = \frac{2J+1}{(2i+1)(2I+1)} \:,
+\end{equation}
+
+\noindent
+and center-of-mass momentum $K_\alpha$ by
+
+\begin{equation}\label{eq:center-of-mass-momentum}
+K_\alpha^2 = (\hbar k_\alpha)^2 = \frac{2mM^2}{(m+M)^2}E \:.
+\end{equation}
+
+\noindent
+Here $E$ is the \textbf{laboratory} kinetic energy of the incident (moving) particle. A derivation of this value for $K_\alpha$ is given in Section II.C.2\ref{}.
+
+The scattering matrix $U$ can be written in terms of matrix $W$ as
+
+\begin{equation}\label{eq:scat-matrix}
+U_{cc'} = \Omega_cW_{cc'}\Omega_{c'} \:,
+\end{equation}
+
+\noindent
+where $\Omega$ is given by
+
+\begin{equation}\label{eq:omega}
+\Omega_c = e^{i(w_c-\phi_c)} \:.
+\end{equation}
+
+\noindent
+Here again, $w_c$ is zero for non-Coulomb channels, and the potential scattering phase shifts for non-Coulomb interactions $\phi_c$ are defined in many references (e.g., \cite{lane_thomas_1958}) and shown in Table \ref{}. The matrix $W$ in Eq. \ref{eq:scat-matrix} is related to the R-matrix (in matrix notation with indices suppressed) via
+
+\begin{equation}\label{eq:W-matrix}
+W = P^{1/2}\left(I-RL\right)^{-1}\left(I-RL^*\right)P^{-1/2} \:.
+\end{equation}
+
+\noindent
+The quantity $I$ in this equation represents the identity matrix, and superscript $*$ indicates a complex conjugate. The form of the R-matrix given in Section IIA.1\ref{} in general Section II.B\ref{} for the versions used in SAMMY. The quantity $L$ in Eq. \ref{eq:W-matrix} is given by
+
+\begin{equation}\label{eq:L-matrix}
+L = (S-B) + iP \:,
+\end{equation}
+
+\noindent
+with $P$ being the penetration factor (penetrability) $S$ the shift factor, and $B$ the arbitrary boundary constant at the channel radius $a_c$. $P$ and $S$ are functions of energy $E$, and also depend on the orbital angular momentum $l$ and the channel radius $a_c$. Formulae for $P$ and $S$ are found in many references (see, for example Eq. (2.9) in \cite{lynn_1968}).
+
+For non-Coulomb interactions, the penetrability and shift factor have the form
+
+\begin{equation}\label{eq:pen-shift-func-of-rho}
+ P\rightarrow P_l(\rho) \qquad \text{and} \qquad S \rightarrow S_l(\rho) \:,
+\end{equation}
+
+\noindent
+where $\rho$ is related to the center-of-mass momentum which in turn is related to the laboratory energy of the incident particle $(E)$. For arbitrary channel $c$ with a particle pair $\alpha$, orbital angular momentum $l$, and channel radius $a_c$, $\rho$ has the form
+
+\begin{equation}\label{eq:de-Broglie-radius}
+ \rho = k_\alpha a_c = \frac{1}{\hbar} \sqrt{\frac{2m_\alpha M_\alpha}{m_\alpha+M_\alpha} \frac{M}{m+M}} \sqrt{(E-\Xi_\alpha)}\: a_c \:,
+\end{equation}
+
+as shown in Section II.C.2\ref{}. Here $\Xi_\alpha$ is the energy threshold for the particle pair $\alpha$, $m_\alpha$ and $M_\alpha$ are the masses of the two particles of particle pair $\alpha$, and $m$ and $M$ are the masses of the incident particle and target nuclide, respectively.
+
+Appropriate formulae for $P$, $S$, and $\phi$ in the non-Coulomb case are shown in Table IIA.1. For two charged particles, formulae for the penetrabilities are given in Section II.C.4\ref{}.
+
+The energy dependence of fission and capture widths is negligible over the energy range of these calculations. Therefore, a penetrability of unity may be used.
+
+% multiline cell: \begin{tabular}{@{}c@{}} line1 \\ line2 \end{tabular}
+\begin{threeparttable}
+\centering
+\caption{Hard-sphere penetrability (penetration factor) $P$, level shift factor $S$, and potential-scattering phase shift $\phi$ for orbital angular momentum $l$, wave number $k$, and channel radius $a_c$, with $\rho=ka_c$} \label{tab:penetrabilities}
+\begin{tabular}{ l c c c }
+ \hline\hline
+ $\mathbf{l}$ & $\mathbf{P_l}$ & $\mathbf{S_l}$ & $\mathbf{\phi_l}$ \cr
+ \hline\hline
+ 0 & $\rho$ & 0 & $\rho$ \cr
+ 1 & $\rho^3/(1+\rho^2)$ & $-1/(1+\rho^2)$ & $\rho-tan^{-1}\rho$ \cr
+ 2 & $\rho^5/(9+3\rho^2+\rho^4)$ & $-(18+3\rho^2)/(9+3\rho^2+\rho^4)$ & $\rho-tan^{-1}\left[3\rho/(3-\rho^2)\right]$ \cr
+ 3 & \begin{tabular}{@{}c@{}}$\rho^7/$ \\ $(225+45\rho^2+6\rho^4+\rho^6)$\end{tabular} & \begin{tabular}{@{}c@{}} $-(675+90\rho^2+6\rho^4)/$ \\ $(225+45\rho^2+6\rho^4+\rho^6)$ \end{tabular} & $\rho-tan^{-1}\left[\rho(15-\rho^2)/(15-6\rho^2)\right]$ \cr
+ \hline
+ $l$ & $\frac{\rho^2P_{l-1}}{(1-S_{l-1})^2+P_{l-1}^2}$ & $\frac{\rho^2(l-S_{l-1})}{(1-S_{l-1})^2+P_{l-1}^2}-l$ & $\phi_{l-1} - tan^{-1}\left(P_{l-1}/(l-S_{l-1})\right)$\tnote{\textdagger} \cr
+
+ \hline\hline
+
+\end{tabular}
+\begin{tablenotes}
+ \item[\textdagger] \footnotesize{The iterative formula for $\phi_l$ could also be defined by $B_l = (B_{l-1}+X_l)/(1-B_{l-1}X_l)$ where $B_l=tan(\rho-\phi_l)$ and $X_l=P_{l-1}/(l-S_l)$}
+\end{tablenotes}
+\end{threeparttable}
+
+
+
+
+
+
+
+
+
diff --git a/sammy/src/blk/CapYCorrections_common.f90 b/sammy/src/blk/Convolution_common.f90
similarity index 69%
rename from sammy/src/blk/CapYCorrections_common.f90
rename to sammy/src/blk/Convolution_common.f90
index e4427788e53c549c8246775c1d2f3eeb8cc8d24a..43c937565f1b94f5c5edd835897b1ccf2d0367fc 100644
--- a/sammy/src/blk/CapYCorrections_common.f90
+++ b/sammy/src/blk/Convolution_common.f90
@@ -3,10 +3,13 @@
! * used to pass the variables contained in C++ classes until
! * the code has been fully modernized.
! **********************************************************/
-module CapYCorrections_common_m
- use CapYieldCorrections_M
+module Convolution_common_m
+ use CapYieldCorrections_M
+ use YieldNormalization_M
IMPLICIT NONE
type(CapYieldCorrections)::capYieldCor
-end module CapYCorrections_common_m
+ type(YieldNormalization)::yNorm
+
+end module Convolution_common_m
diff --git a/sammy/src/blk/Observable_common.f90 b/sammy/src/blk/Observable_common.f90
new file mode 100644
index 0000000000000000000000000000000000000000..7fd5eb21f37dcd6ff1950d9ec763ddd3e1becfd7
--- /dev/null
+++ b/sammy/src/blk/Observable_common.f90
@@ -0,0 +1,12 @@
+! /**********************************************************
+! * This common block style module is a temporary solution
+! * used to pass the variables contained in C++ classes until
+! * the code has been fully modernized.
+! **********************************************************/
+module Observable_common_m
+ use Yield_M
+ IMPLICIT NONE
+
+ type(Yield)::yld
+
+end module Observable_common_m
diff --git a/sammy/src/convolution/CMakeLists.txt b/sammy/src/convolution/CMakeLists.txt
index 81c874b5866af30d9dca4b5ca8b89f7ac1056214..f576ac814532c1f91fe6e8075fc7092eaacf1259 100644
--- a/sammy/src/convolution/CMakeLists.txt
+++ b/sammy/src/convolution/CMakeLists.txt
@@ -11,8 +11,10 @@ INCLUDE_DIRECTORIES(${CMAKE_CURRENT_BINARY_DIR})
# Set headers
SET(HEADERS
CapYieldCorrections.h
-
interface/cpp/CapYieldCorrectionsInterface.h
+
+ YieldNormalization.h
+ interface/cpp/YieldNormalizationInterface.h
)
APPEND_SET(CONVOLUTION_SOURCES
@@ -21,6 +23,11 @@ APPEND_SET(CONVOLUTION_SOURCES
interface/cpp/CapYieldCorrectionsInterface.cpp
interface/fortran/CapYieldCorrections_I.f90
interface/fortran/CapYieldCorrections_M.f90
+
+ YieldNormalization.cpp
+ interface/cpp/YieldNormalizationInterface.cpp
+ interface/fortran/YieldNormalization_I.f90
+ interface/fortran/YieldNormalization_M.f90
)
diff --git a/sammy/src/convolution/CapYieldCorrections.cpp b/sammy/src/convolution/CapYieldCorrections.cpp
index ace8fdc6c18a716b1f3b10e52be46f6ef7c22eab..8c8f736accb253b7b47ea6777f9890651b160951 100644
--- a/sammy/src/convolution/CapYieldCorrections.cpp
+++ b/sammy/src/convolution/CapYieldCorrections.cpp
@@ -86,6 +86,9 @@ namespace sammy{
if( N_sigma <= 0.0 ){
std::ostringstream err;
err << "(sample thickness [at/b]) * (tot. x.s. [b]) is less than zero" << std::endl;
+ err << "--------------------------------" << std::endl;
+ err << "N * sigma_tot = " << N_sigma << std::endl;
+ err << "--------------------------------" << std::endl;
throw std::runtime_error(err.str());
}
@@ -94,10 +97,11 @@ namespace sammy{
std::ostringstream err;
err << "Under or overflow may occur in exponent. Check sample thickness" << std::endl;
err << " and total cross section values." << std::endl;
+ err << "--------------------------------" << std::endl;
+ err << "N * sigma_tot = " << N_sigma << std::endl;
+ err << "--------------------------------" << std::endl;
throw std::runtime_error(err.str());
}
-
-
}
void CapYieldCorrections::check_attenFactor(const double &attenFactor)
@@ -106,7 +110,9 @@ namespace sammy{
if( attenFactor <= 0.0 ){
std::ostringstream err;
err << "Gamma attenuation factor C1 (from documentation) is less than zero" << std::endl;
+ err << "--------------------------------" << std::endl;
err << "C1 = " << attenFactor << std::endl;
+ err << "--------------------------------" << std::endl;
throw std::runtime_error(err.str());
}
@@ -115,7 +121,9 @@ namespace sammy{
std::ostringstream err;
err << "Gamma attenuation factor C1 (from documentation) is too large, and" << std::endl;
err << " may cause under/overflow issues." << std::endl;
+ err << "--------------------------------" << std::endl;
err << "C1 = " << attenFactor << std::endl;
+ err << "--------------------------------" << std::endl;
throw std::runtime_error(err.str());
}
}
@@ -126,9 +134,10 @@ namespace sammy{
std::ostringstream err;
err << "The gamma attenuation correction to the yield is less than or" << std::endl;
err << "equal to zero or greater than 1.0" << std::endl;
+ err << "--------------------------------" << std::endl;
err << "k_a = " << gammaCor << std::endl;
+ err << "--------------------------------" << std::endl;
throw std::runtime_error(err.str());
}
}
-
}
\ No newline at end of file
diff --git a/sammy/src/convolution/YieldNormalization.cpp b/sammy/src/convolution/YieldNormalization.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..1c0ed3720e640818047cbed582779e2e1ace6b10
--- /dev/null
+++ b/sammy/src/convolution/YieldNormalization.cpp
@@ -0,0 +1,142 @@
+
+#include "YieldNormalization.h"
+
+#include <iostream>
+namespace sammy{
+
+ YieldNormalization::YieldNormalization(){}
+
+ void YieldNormalization::norm(sammy::Yield & yld, sammy::SampleDimensions & sd,
+ sammy::MultipleScattering &multScat, double & sigTot)
+ {
+ // How many scatters considered before neutron absorption?
+ // (We change the behavior of normTypes == 1 and 2 for different scatters)
+ int NS = multScat.getNumScatters();
+
+ if( yld.normTypeIsSelfShielded() )
+ {
+ selfShielded(yld,sd);
+ }
+ else if( yld.normTypeIsDivideByN() )
+ {
+ divideByN(yld, sd, NS);
+ }
+ else if( yld.normTypeIsTimesSigTot() )
+ {
+ timesSigmaTot(yld, sd, sigTot, NS);
+ }
+ else
+ {
+ std::ostringstream err;
+ err << "------------------------------------------" << std::endl;
+ err << "The yield normalization type has not been set, check input file." << std::endl;
+ err << "------------------------------------------" << std::endl;
+ throw std::runtime_error(err.str());
+ }
+ }
+
+ void YieldNormalization::selfShielded(sammy::Yield &yld, sammy::SampleDimensions &sd)
+ {
+ yld.setY0( yld.getY0() * sd.getArealDensity() );
+ }
+
+ void YieldNormalization::divideByN(sammy::Yield &yld, sammy::SampleDimensions &sd, int numScatters)
+ {
+ if( numScatters == 0 )
+ {
+ // do nothing
+ }
+
+ if( numScatters == 1 or numScatters == 2 )
+ {
+ // Yield
+ yld.setY1( yld.getY1() / sd.getArealDensity() );
+
+ // if numParam == 0, we don't need to calculate
+ for( int i=0; i<yld.getNumTotalParam(); i++ )
+ {
+ // Derivative of yield
+ yld.setDerivY1(i, yld.getDerivY1(i) / sd.getArealDensity() );
+ }
+ }
+
+ // Really if the number of scatters is >= 2
+ if( numScatters == 2 )
+ {
+ // Yield
+ yld.setY2( yld.getY2() / sd.getArealDensity() );
+
+ // if numParam == 0, we don't need to calculate
+ for( int i=0; i<yld.getNumTotalParam(); i++ )
+ {
+ // Derivative of yield
+ yld.setDerivY2(i, yld.getDerivY2(i) / sd.getArealDensity() );
+ }
+ }
+ }
+
+ void YieldNormalization::timesSigmaTot(sammy::Yield &yld, sammy::SampleDimensions &sd,
+ double & sigTot, int numScatters)
+ {
+ // for all cases
+ yld.setY0( yld.getY0() * sigTot * sd.getArealDensity() );
+
+ if( numScatters == 0 )
+ {
+ // nothing extra
+ }
+
+ if( numScatters == 1 or numScatters == 2 )
+ {
+ yld.setY1( yld.getY1() * sigTot );
+
+ // if numParam == 0, we don't need to calculate
+ for( int i=0; i<yld.getNumTotalParam(); i++ )
+ {
+ yld.setDerivY1(i, yld.getDerivY1(i) * sigTot );
+ }
+ }
+
+ if( numScatters == 2 )
+ {
+ yld.setY2( yld.getY2() * sigTot );
+
+ // if numParam == 0, we don't need to calculate
+ for( int i=0; i<yld.getNumTotalParam(); i++ )
+ {
+ yld.setDerivY2(i, yld.getDerivY2(i) * sigTot );
+ }
+ }
+ }
+
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/sammy/src/convolution/YieldNormalization.h b/sammy/src/convolution/YieldNormalization.h
new file mode 100644
index 0000000000000000000000000000000000000000..b592d3f834393ca635e86ba130eae6a9b0ba24a1
--- /dev/null
+++ b/sammy/src/convolution/YieldNormalization.h
@@ -0,0 +1,60 @@
+#ifndef YIELDNORMALIZATION_H
+#define YIELDNORMALIZATION_H
+
+#include "../salmon/Yield.h"
+#include "../salmon/ExperimentalParameters.h"
+#include "../salmon/MultipleScattering.h"
+
+namespace sammy{
+
+ /** @brief Normalization of the observed yield
+
+ This class applies normalization to yields (capture, fission, scattering...)
+ @date February 2020
+ */
+ class YieldNormalization
+ {
+ public:
+ YieldNormalization();
+ YieldNormalization(const YieldNormalization& orig);
+ virtual ~YieldNormalization(){};
+
+ /**
+ * Master normalization function
+ */
+ void norm(sammy::Yield & yld, sammy::SampleDimensions & sd,
+ sammy::MultipleScattering &multScat, double & sigTot );
+
+ /**
+ * Normalize as self-shielded yield \f$ Y \f$
+ * \f[ = \frac{\sigma_x}{\sigma_t} \left(1-e^{-n\sigma_t}\right) \f]
+ */
+ void selfShielded(sammy::Yield & yld, sammy::SampleDimensions &sd);
+
+ /**
+ * Normalize as self-shielded yield divided by sample thickness \f$ Y/n \f$ which, in the
+ * thin sample approximation, is \f$ \approx \sigma_x \f$
+ * \f[ = \frac{\sigma_\gamma}{\sigma_t n} \left(1-e^{-n\sigma_t}\right) \f]
+ */
+ void divideByN(sammy::Yield & yld, sammy::SampleDimensions &sd, int numScatters);
+
+ /**
+ * Normalize as self-shielded yield \f$ Y \f$ multiplied by total cross section \f$ \sigma_t \f$
+ * \f[ = \sigma_x \left(1-e^{-n\sigma_t}\right) \f]
+ */
+ void timesSigmaTot(sammy::Yield & yld, sammy::SampleDimensions &sd,
+ double & sigTot, int numScatters);
+
+ /**
+ * Check that the number of varied parameters is greater than zero if
+ * Bayes' equations are to be solved
+ * @param yld A sammy::Yield object
+ */
+ //void checkNumParam(sammy::Yield & yld);
+
+ };
+
+}
+
+
+#endif // YIELDNORMALIZATION_H
\ No newline at end of file
diff --git a/sammy/src/convolution/cmake/Dependencies.cmake b/sammy/src/convolution/cmake/Dependencies.cmake
index 13244a48ada913fc64d9cc8d62923c6e8cffd7a1..1cbd5653294e182b98d8c4f58fced72846583024 100644
--- a/sammy/src/convolution/cmake/Dependencies.cmake
+++ b/sammy/src/convolution/cmake/Dependencies.cmake
@@ -1,4 +1,4 @@
-SET(LIB_REQUIRED_DEP_PACKAGES)
+SET(LIB_REQUIRED_DEP_PACKAGES salmon)
SET(LIB_OPTIONAL_DEP_PACKAGES)
SET(TEST_REQUIRED_DEP_PACKAGES)
SET(TEST_OPTIONAL_DEP_PACKAGES)
diff --git a/sammy/src/convolution/interface/cix/YieldNormalization.cpp2f.xml b/sammy/src/convolution/interface/cix/YieldNormalization.cpp2f.xml
new file mode 100644
index 0000000000000000000000000000000000000000..c3fc0a03e01190ee4d75e013ee37ed894da09b11
--- /dev/null
+++ b/sammy/src/convolution/interface/cix/YieldNormalization.cpp2f.xml
@@ -0,0 +1,15 @@
+<generate name="YieldNormalization">
+ <include_relative name="../../YieldNormalization.h"/>
+ <using_namespace name="sammy"/>
+
+ <class name="YieldNormalization">
+
+ <method name="norm">
+ <param name="yld" type="Yield"/>
+ <param name="sampDim" type="SampleDimensions"/>
+ <param name="multscat" type="MultipleScattering"/>
+ <param name="sigTot" type="double"/>
+ </method>
+
+ </class>
+</generate>
\ No newline at end of file
diff --git a/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.cpp b/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..908f0efc931ebd738660fc35926c4b1023ceca8a
--- /dev/null
+++ b/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.cpp
@@ -0,0 +1,26 @@
+/*!
+* This file has been dynamically generated by Class Interface Xml (CIX)
+* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+* If changes need to occur, modify the appropriate CIX xml file
+* Date Generated: Mon Mar 02 10:18:42 EST 2020
+* If any issues are experiences with this generated file that cannot be fixed
+* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+*/
+#include <string.h>
+#include "YieldNormalizationInterface.h"
+using namespace sammy;
+void YieldNormalization_norm(void * YieldNormalization_ptr,Yield * yld,SampleDimensions * sampDim,MultipleScattering * multscat,double * sigTot)
+{
+ ((YieldNormalization*)YieldNormalization_ptr)->norm(*yld,*sampDim,*multscat,*sigTot);
+}
+
+void* YieldNormalization_initialize()
+{
+ return new YieldNormalization();
+}
+
+void YieldNormalization_destroy(void * YieldNormalization_ptr)
+{
+ delete (YieldNormalization*)YieldNormalization_ptr;
+}
+
diff --git a/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.h b/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.h
new file mode 100644
index 0000000000000000000000000000000000000000..36db444e93f859d5cf28d5ed024353f2df23a8f3
--- /dev/null
+++ b/sammy/src/convolution/interface/cpp/YieldNormalizationInterface.h
@@ -0,0 +1,22 @@
+/*!
+* This file has been dynamically generated by Class Interface Xml (CIX)
+* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+* If changes need to occur, modify the appropriate CIX xml file
+* Date Generated: Mon Mar 02 10:18:42 EST 2020
+* If any issues are experiences with this generated file that cannot be fixed
+* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+*/
+#ifndef YIELDNORMALIZATIONINTERFACE_H
+#define YIELDNORMALIZATIONINTERFACE_H
+#include "../../YieldNormalization.h"
+using namespace sammy;
+#ifdef __cplusplus
+extern "C" {
+#endif
+void YieldNormalization_norm(void * YieldNormalization_ptr,Yield * yld,SampleDimensions * sampDim,MultipleScattering * multscat,double * sigTot);
+void* YieldNormalization_initialize();
+void YieldNormalization_destroy(void * YieldNormalization_ptr);
+#ifdef __cplusplus
+}
+#endif
+#endif /* YIELDNORMALIZATIONINTERFACE_H */
diff --git a/sammy/src/convolution/interface/fortran/YieldNormalization_I.f90 b/sammy/src/convolution/interface/fortran/YieldNormalization_I.f90
new file mode 100644
index 0000000000000000000000000000000000000000..4610ddb66eeb4b99807e03ddc3dcb6dd311fc26e
--- /dev/null
+++ b/sammy/src/convolution/interface/fortran/YieldNormalization_I.f90
@@ -0,0 +1,32 @@
+!>
+!! This file has been dynamically generated by Class Interface Xml (CIX)
+!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+!! If changes need to occur, modify the appropriate CIX xml file
+!! Date Generated: Mon Mar 02 10:18:42 EST 2020
+!! If any issues are experiences with this generated file that cannot be fixed
+!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+!!/
+module YieldNormalization_I
+use, intrinsic :: ISO_C_BINDING
+interface
+subroutine f_YieldNormalization_norm(YieldNormalization_ptr, yld,sampDim,multscat,sigTot ) BIND(C,name="YieldNormalization_norm")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: YieldNormalization_ptr;
+ type(C_PTR), value :: yld;
+ type(C_PTR), value :: sampDim;
+ type(C_PTR), value :: multscat;
+ real(C_DOUBLE) :: sigTot;
+end subroutine
+type(C_PTR) function f_YieldNormalization_initialize( )BIND(C,name="YieldNormalization_initialize")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR) :: YieldNormalization_ptr;
+end function
+subroutine f_YieldNormalization_destroy(this) BIND(C,name="YieldNormalization_destroy")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: this;
+end subroutine
+end interface
+end module YieldNormalization_I
diff --git a/sammy/src/convolution/interface/fortran/YieldNormalization_M.f90 b/sammy/src/convolution/interface/fortran/YieldNormalization_M.f90
new file mode 100644
index 0000000000000000000000000000000000000000..0473ffa31ee6d258b03b5b14fc71de85fc8d902d
--- /dev/null
+++ b/sammy/src/convolution/interface/fortran/YieldNormalization_M.f90
@@ -0,0 +1,43 @@
+!>
+!! This file has been dynamically generated by Class Interface Xml (CIX)
+!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+!! If changes need to occur, modify the appropriate CIX xml file
+!! Date Generated: Mon Mar 02 10:18:42 EST 2020
+!! If any issues are experiences with this generated file that cannot be fixed
+!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+!!/
+module YieldNormalization_M
+use, intrinsic :: ISO_C_BINDING
+use YieldNormalization_I
+use MultipleScattering_M
+use ExperimentalParameters_M
+use Yield_M
+type YieldNormalization
+ type(C_PTR) :: instance_ptr=C_NULL_PTR
+ contains
+ procedure, pass(this) :: norm => YieldNormalization_norm
+ procedure, pass(this) :: initialize => YieldNormalization_initialize
+ procedure, pass(this) :: destroy => YieldNormalization_destroy
+end type YieldNormalization
+contains
+subroutine YieldNormalization_norm(this, yld, sampDim, multscat, sigTot)
+ implicit none
+ class(YieldNormalization)::this
+ class(Yield)::yld
+ class(SampleDimensions)::sampDim
+ class(MultipleScattering)::multscat
+ real(C_DOUBLE)::sigTot
+ call f_YieldNormalization_norm(this%instance_ptr, yld%instance_ptr,sampDim%instance_ptr,multscat%instance_ptr,sigTot)
+end subroutine
+subroutine YieldNormalization_initialize(this)
+ implicit none
+ class(YieldNormalization) :: this
+ this%instance_ptr = f_YieldNormalization_initialize()
+end subroutine
+subroutine YieldNormalization_destroy(this)
+ implicit none
+ class(YieldNormalization) :: this
+ call f_YieldNormalization_destroy(this%instance_ptr)
+ this%instance_ptr = C_NULL_PTR
+end subroutine
+end module YieldNormalization_M
diff --git a/sammy/src/convolution/tests/CMakeLists.txt b/sammy/src/convolution/tests/CMakeLists.txt
index 50399bf915562fe3fb2657ece7f75a5a29fc9194..1fb53c9229e8f3e29f4d3801a12fa1a9e2bdf094 100644
--- a/sammy/src/convolution/tests/CMakeLists.txt
+++ b/sammy/src/convolution/tests/CMakeLists.txt
@@ -1,3 +1,4 @@
include(NemesisTest)
ADD_NEMESIS_TEST(CapYieldCorrectionsTest.cpp TIMEOUT 60000 NP 1 )
+ADD_NEMESIS_TEST(YieldNormalizationTest.cpp TIMEOUT 60000 NP 1 )
diff --git a/sammy/src/convolution/tests/YieldNormalizationTest.cpp b/sammy/src/convolution/tests/YieldNormalizationTest.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..6dd8ee1a9a1eb36a28483178420c1827d84fe42d
--- /dev/null
+++ b/sammy/src/convolution/tests/YieldNormalizationTest.cpp
@@ -0,0 +1,353 @@
+#include "Nemesis/gtest/Gtest_Functions.hh"
+#include "Nemesis/gtest/nemesis_gtest.hh"
+#include "../YieldNormalization.h"
+#include <exception>
+#include <stdexcept>
+
+using namespace sammy;
+
+/**************************************************************************************************
+* TODO: need to include tests for the Observable class members of yield from Observable<Yield>
+**************************************************************************************************/
+
+/**************************************************************************************************
+* Test norm function: multiply by tot x.s., scattering = 0
+**************************************************************************************************/
+TEST( YieldNormalization,multTot_scat0 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeTimesSigTot();
+ ms.setNumScatters(0);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setDerivY0(0,-0.01);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * sigTot * 0.01;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+
+}
+
+/**************************************************************************************************
+* Test norm function: multiply by tot x.s., scattering = 1
+**************************************************************************************************/
+TEST( YieldNormalization,multTot_scat1 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeTimesSigTot();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(1);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * sigTot * 0.01;
+ double gold_Y1 = 0.2 * sigTot;
+ double gold_DerivY1 = -0.02 * sigTot;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: multiply by tot x.s., scattering = 2
+**************************************************************************************************/
+TEST( YieldNormalization,multTot_scat2 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeTimesSigTot();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(2);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ yld.setY2(0.3);
+ yld.setDerivY2(0,-0.03);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * sigTot * 0.01;
+ double gold_Y1 = 0.2 * sigTot;
+ double gold_Y2 = 0.3 * sigTot;
+
+ double gold_DerivY1 = -0.02 * sigTot;
+ double gold_DerivY2 = -0.03 * sigTot;
+
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_Y2, yld.getY2());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+ ASSERT_EQ(gold_DerivY2,yld.getDerivY2(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: divide by thick, scattering = 0
+**************************************************************************************************/
+TEST( YieldNormalization,divideByN_scat0 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeDivideByN();
+ ms.setNumScatters(0);
+
+ yld.setY0(0.1);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+
+}
+
+/**************************************************************************************************
+* Test norm function: divide by thick, scattering = 1
+**************************************************************************************************/
+TEST( YieldNormalization,divideByN_scat1 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeDivideByN();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(1);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1;
+ double gold_Y1 = 0.2 / 0.01;
+ double gold_DerivY1 = -0.02 / 0.01;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: divide by thick, scattering = 2
+**************************************************************************************************/
+TEST( YieldNormalization,divideByN_scat2 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeDivideByN();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(2);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ yld.setY2(0.3);
+ yld.setDerivY2(0,-0.03);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1;
+ double gold_Y1 = 0.2 / 0.01;
+ double gold_DerivY1 = -0.02 / 0.01;
+ double gold_Y2 = 0.3 / 0.01;
+ double gold_DerivY2 = -0.03 / 0.01;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+ ASSERT_EQ(gold_Y2, yld.getY2());
+ ASSERT_EQ(gold_DerivY2,yld.getDerivY2(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: self-shielded norm, scattering = 0
+**************************************************************************************************/
+TEST( YieldNormalization,selfShielded_scat0 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeSelfShielded();
+ ms.setNumScatters(0);
+
+ yld.setY0(0.1);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * 0.01;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+
+}
+
+/**************************************************************************************************
+* Test norm function: self-shielded norm, scattering = 1
+**************************************************************************************************/
+TEST( YieldNormalization,selfShielded_scat1 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeSelfShielded();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(1);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * 0.01;
+ double gold_Y1 = 0.2;
+ double gold_DerivY1 = -0.02;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: self-shielded norm, scattering = 2
+**************************************************************************************************/
+TEST( YieldNormalization,selfShielded_scat2 ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ YieldNormalization yNorm;
+
+ // These change test-to-test
+ yld.setNormTypeSelfShielded();
+ yld.setCalcDerivs(true);
+ ms.setNumScatters(2);
+ yld.setNumTotalParam(1);
+
+ yld.setY0(0.1);
+ yld.setY1(0.2);
+ yld.setDerivY1(0,-0.02);
+ yld.setY2(0.3);
+ yld.setDerivY2(0,-0.03);
+ double sigTot = 3.0;
+ sd.setArealDensity(0.01);
+
+ yNorm.norm(yld, sd, ms, sigTot);
+
+ double gold_Y0 = 0.1 * 0.01;
+ double gold_Y1 = 0.2;
+ double gold_DerivY1 = -0.02;
+ double gold_Y2 = 0.3;
+ double gold_DerivY2 = -0.03;
+
+ ASSERT_EQ(gold_Y0, yld.getY0());
+ ASSERT_EQ(gold_Y1, yld.getY1());
+ ASSERT_EQ(gold_DerivY1,yld.getDerivY1(0));
+ ASSERT_EQ(gold_Y2, yld.getY2());
+ ASSERT_EQ(gold_DerivY2,yld.getDerivY2(0));
+
+}
+
+/**************************************************************************************************
+* Test norm function: normalization not set, throw the user an error
+**************************************************************************************************/
+TEST( YieldNormalization,normNotSet ){
+
+ Yield yld;
+ SampleDimensions sd;
+ MultipleScattering ms;
+
+ double sigTot = 3.0;
+
+ YieldNormalization yNorm;
+
+ // Expect to throw error
+ int gold_result = 0;
+ int result;
+
+ try{
+ yNorm.norm(yld, sd, ms, sigTot);
+ result = 1;
+ }
+ catch(...){
+ result = 0;
+ }
+
+
+ ASSERT_EQ(gold_result,result);
+}
diff --git a/sammy/src/endf/CovarianceData.h b/sammy/src/endf/CovarianceData.h
index 9c1031ea55a6de0343987128e419066ad5f517ac..602ae91a6c27865e52e57514569a97ffc38fba79 100644
--- a/sammy/src/endf/CovarianceData.h
+++ b/sammy/src/endf/CovarianceData.h
@@ -173,7 +173,7 @@ namespace sammy{
/**
- * Get the number of total parameters
+ * Get the number of total parameters: The pup'd and non-pup'd
*/
int getNumTotalParam() const{
return (int)paramIndex.size();
diff --git a/sammy/src/inp/minp04.F b/sammy/src/inp/minp04.F
index e58f420c31dbc7a860deb338bd731fdefec2c2aa..2e68ef8e7595d33ed16525aae7ab1a629d9b9bb3 100644
--- a/sammy/src/inp/minp04.F
+++ b/sammy/src/inp/minp04.F
@@ -14,6 +14,7 @@ C
use lbro_common_m
use misccc_common_m
use ntyp_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
C
C
@@ -81,25 +82,26 @@ C
C
IF (Kdebug.EQ.0) Debug = .False.
IF (Kdebug.EQ.1) Debug = .True.
- IF (Kyield.EQ.-1 .AND. (Kssmsc.GT.0 .OR. Kssmsc.EQ.-2)) THEN
- WRITE (6,10000)
- WRITE (21,10000)
-10000 FORMAT (/, ' ### CAUTION -- Normalization method has not been
- * specified in INPut file. ###', /,
- * ' ### See page 56r of the SAMMY users manual for a des
- *cription ###', /,
- * ' ### of the various options. Please add the appropria
- *te command ###', /,
- * ' ### to your INPut file and resubmit. Possibilities x
- *include ###', /,
- * ' ### NORMALIZE AS CROSS SECTION or CROSS SECTION x
- * ###', /,
- * ' ### NORMALIZE AS YIELD RATHER or YIELD x
- * ###', /,
- * ' ### NORMALIZE AS (1-E)SIGMA x
- * ###')
- STOP
- END IF
+C IF (yld%getNormType().EQ.-1 .AND.
+C * (Kssmsc.GT.0 .OR. Kssmsc.EQ.-2)) THEN
+C WRITE (6,10000)
+C WRITE (21,10000)
+C 10000 FORMAT (/, ' ### CAUTION -- Normalization method has not been
+C * specified in INPut file. ###', /,
+C * ' ### See page 56r of the SAMMY users manual for a des
+C *cription ###', /,
+C * ' ### of the various options. Please add the appropria
+C *te command ###', /,
+C * ' ### to your INPut file and resubmit. Possibilities x
+C *include ###', /,
+C * ' ### NORMALIZE AS CROSS SECTION or CROSS SECTION x
+C * ###', /,
+C * ' ### NORMALIZE AS YIELD RATHER or YIELD x
+C * ###', /,
+C * ' ### NORMALIZE AS (1-E)SIGMA x
+C * ###')
+C STOP
+C END IF
C
IF (Mretro.EQ.1) THEN
IF (Nretro.EQ.1) THEN
@@ -189,6 +191,7 @@ C
use broad_common_m
use namfil_common_m
use lbro_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
C
#ifdef COMPILE_WITH_SAMINT
@@ -533,11 +536,12 @@ C
C *** Do not make plot file of multiple-scattering pieces (-268,269)
Kwssms = 0
C
-C *** Normalize as cross section rather than yield
+C *** Normalize as cross section rather than yield
+C *** Kyield=yld%getNormType()
C *** for Kyield=1: (270,271) cross section
C *** for Kyield=0: (272,273) yield
C *** for Kyield=2: (274 ) 1-e
- Kyield = -1
+ !call yld%setNormType(-1)
C
C *** Do not prepare ssm comparison file (-275)
C *** (i.e., do not print multiple-scattering corrections)
diff --git a/sammy/src/inp/minp05.F b/sammy/src/inp/minp05.F
index 2c6bd23e0f59f3bca7e10ed3b0ce3a426322ef01..6292686e528b14afb54618ed572ed615b444f8c3 100644
--- a/sammy/src/inp/minp05.F
+++ b/sammy/src/inp/minp05.F
@@ -15,6 +15,7 @@ C
use namfil_common_m
use lbro_common_m
use misccc_common_m
+ use Observable_common_m
use mdf5_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
@@ -416,11 +417,11 @@ C
ELSE IF (Number.EQ.268 .OR. Number.EQ.269) THEN
Kwssms = 1
ELSE IF (Number.EQ.270 .OR. Number.EQ.271) THEN
- Kyield = 1
+ call yld%setNormTypeDivideByN()
ELSE IF (Number.EQ.272 .OR. Number.EQ.273) THEN
- Kyield = 0
+ call yld%setNormTypeSelfShielded()
ELSE IF (Number.EQ.274) THEN
- Kyield = 2
+ call yld%setNormTypeTimesSigTot()
ELSE IF (Number.EQ.275) THEN
Kssmpr = 1
ELSE IF (Number.EQ.276 .OR. Number.eq.277) THEN
diff --git a/sammy/src/mso/mmso5.f b/sammy/src/mso/mmso5.f
index 792fd5a5fa6ebee8f9ddcc1809b7cd21a6ee5d56..23fff9a61676887840f8b075c09c778ae6484e30 100644
--- a/sammy/src/mso/mmso5.f
+++ b/sammy/src/mso/mmso5.f
@@ -487,13 +487,14 @@ C
use fixedi_m
use ifwrit_m
use fixedr_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
DIMENSION Capsig(*), Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*),
* Dy2(*), Sigxxx(*), Dasigx(*), Dbsigx(*)
DATA Zero /0.0d0/, Half /0.5d0/, One /1.0d0/
C
C *** Default normalization is such that Y0 = (1-e)Capture /(n Total)
-C *** Which corresponds to Kyield=1
+C *** Which corresponds to yld%normTypeIsDivideByN()==true
C
Cap = Capsig(Nn)
IF (Kssmsc.NE.-2) THEN
@@ -522,15 +523,15 @@ C
Sigxxx(1) = Y0
IF (Kssmsc.NE.-2) Sigxxx(1) = Sigxxx(1) + Y1
IF (Kssdbl.EQ. 1) Sigxxx(1) = Sigxxx(1) + Y2
- IF (Kyield.EQ. 0) Sigxxx(1) = Sigxxx(1)*Dthick
- IF (Kyield.EQ. 2) Sigxxx(1) = Sigxxx(1)*Total*Dthick
+ IF (yld%normTypeIsSelfShielded() ) Sigxxx(1) = Sigxxx(1)*Dthick
+ IF (yld%normTypeIsTimesSigTot()) Sigxxx(1) =Sigxxx(1)*Total*Dthick
C
IF (Kssmpr.EQ.1) THEN
- IF (Kyield.EQ.0) THEN
+ IF ( yld%normTypeIsSelfShielded() ) THEN
Y0 = Y0 * Dthick
Y1 = Y1 * Dthick
Y2 = Y2 * Dthick
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF ( yld%normTypeIsTimesSigTot() ) THEN
Y0 = Y0 * Dthick * Total
Y1 = Y1 * Dthick * Total
Y2 = Y2 * Dthick * Total
@@ -551,9 +552,9 @@ C
IF (Kvthck.GT.0) THEN
C *** Add deriv of Y0 wrt Thickness; modify deriv of Y1 & Y2 wrt Thick
Kv = Kvthck - Nvadif - Ndasig
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Dbsigx(Kv) = Exp1*Cap + Dbsigx(Kv)
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Sn = Dthick*Total
IF (Sn.LT.0.3d0) THEN
Abc = Abcexp (-Sn, A, B, C, D, Ijklmn)
@@ -564,14 +565,14 @@ C *** Add deriv of Y0 wrt Thickness; modify deriv of Y1 & Y2 wrt Thick
Dbsigx(Kv) = F + Dbsigx(Kv)
IF (Kssmsc.NE.-2) Dbsigx(Kv) = Dbsigx(Kv) - Y1/Dthick
IF (Kssdbl.EQ. 1) Dbsigx(Kv) = Dbsigx(Kv) - Y2/Dthick
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Dbsigx(Kv)= Exp1*Cap*Total + Dbsigx(Kv)
END IF
END IF
C
IF (Sensin.NE.Zero .AND. Ksensc.GT.0) THEN
- IF (Kyield.EQ.0) Asensn = Asensn * Dthick
- IF (Kyield.EQ.2) Asensn = Asensn * Dthick * Total
+ IF (yld%normTypeIsSelfShielded()) Asensn = Asensn * Dthick
+ IF (yld%normTypeIsTimesSigTot()) Asensn = Asensn *Dthick *Total
C *** Add deriv of Y0 wrt neutron sensitivity multiplier
Kv = Ksensc - Nvadif - Ndasig
Dbsigx(Kv) = Asensn
@@ -586,6 +587,7 @@ C
* Dy2, Y0, Y1, Y2, Dthick, Total, Exp1, Dddyy0, Nn, Nd, Na)
use fixedi_m
use ifwrit_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
DIMENSION Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*), Dy2(*),ds(*)
DATA One /1.0d0/
@@ -594,40 +596,40 @@ C
Iipar = Ipar + Na
C
C *** First, the derivative of Y0 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) +
* (Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total)
* + Dtotsi(Iipar,Nn)*Exp1*Cap/Total ) * Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) +
* Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total)
* + Dtotsi(Iipar,Nn)*Exp1*Cap/Total
IF (Iipar.EQ.1) Dddyy0 = Dddyy0 +
* Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total)
* + Dtotsi(Iipar,Nn)*Exp1*Cap/Total
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + (One-Exp1)*Dcapsi(Iipar,Nn)
* + Dtotsi(Iipar,Nn)*Exp1*Cap*Dthick
END IF
C
IF (Kssmsc.NE.-2) THEN
C *** Next, the derivative of Y1 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)*Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)*Total*Dthick +
* Y1*Dtotsi(Iipar,Nn)*Dthick
END IF
C
IF (Kssdbl.EQ.1) THEN
C *** Finally, the derivative of Y2 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar)*Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar)*Total*Dthick
* + Y2*Dtotsi(Iipar,Nn)*Dthick
END IF
diff --git a/sammy/src/par/mpar0.f b/sammy/src/par/mpar0.f
index afaa8d4456cb71dc59c9a9bf141770fc1766029c..3650a3b7ec7f3ce5e265ca8b5a84081b7ba39f84 100644
--- a/sammy/src/par/mpar0.f
+++ b/sammy/src/par/mpar0.f
@@ -427,7 +427,9 @@ C
use fixedr_m
use ExpPars_common_m
use EndfData_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
+ character(len=80)::norm_string
DATA Zero /0.0d0/
real(kind=8) Ybeam
@@ -448,17 +450,15 @@ C
10700 FORMAT ('sigma.dat')
WRITE (12,10800) Lllmax-1, Numiso
10800 FORMAT (2I3)
- IF (Kyield.EQ.0) THEN
- WRITE (12,10900)
-10900 FORMAT ('Yield = Normalization type')
- ELSE IF (Kyield.EQ.1) THEN
- WRITE (12,11000)
-11000 FORMAT ('Cross section = Normalization type')
- ELSE IF (Kyield.EQ.2) THEN
- WRITE (12,11100)
-11100 FORMAT ('Other (1-e)sigma = Normalization type')
+ call yld%getNormString(norm_string)
+ IF ( yld%normTypeIsSelfShielded() ) THEN
+ WRITE (12,*) norm_string
+ ELSE IF ( yld%normTypeIsDivideByN() ) THEN
+ WRITE (12,*) norm_string
+ ELSE IF ( yld%normTypeIsTimesSigTot() ) THEN
+ WRITE (12,*) norm_string
ELSE
- STOP '[STOP in mpar0.f, Kyield is unknown]'
+ STOP '[STOP in mpar0.f, Kyield/( normaliz. type) is unknown]'
END IF
IF (Numiso.GT.0) THEN
WRITE (12,11200) (resParData%getMassForIsotope(I),I=1,Numiso)
diff --git a/sammy/src/ref/mcon.F b/sammy/src/ref/mcon.F
index 2708283934e6627d0a1d8334fad5f5881170646a..7a9852c995fb6249aadcaeac7b4b68000d7bf782 100644
--- a/sammy/src/ref/mcon.F
+++ b/sammy/src/ref/mcon.F
@@ -19,6 +19,7 @@ C These are here to replace the common block
use zzzzz_common_m
use EndfData_common_m
use Samdat_0_M
+ use Observable_common_m
C this is not a common block, this contains functions
use mold4_m
IMPLICIT DOUBLE PRECISION (A-h,o-z)
@@ -31,6 +32,7 @@ CC ****************************************************
C
real(kind=8)::A(-1*(SAMMY_ARRAY_SIZE):(SAMMY_ARRAY_SIZE))
+ call yld%initialize()
call resParData%initialize()
call covData%initialize()
call covData%makeNewCovariance()
diff --git a/sammy/src/salmon/CMakeLists.txt b/sammy/src/salmon/CMakeLists.txt
index f80412a95af0eb7ae4fcb3058037ce7075ac21dd..c26baf928b07c490db21250421a712baf0cb412b 100644
--- a/sammy/src/salmon/CMakeLists.txt
+++ b/sammy/src/salmon/CMakeLists.txt
@@ -1,5 +1,8 @@
INCLUDE(SammyPackageSetup)
+# --------------------------------------
+# SAmmy Legacy Memory Obsolete Now
+# --------------------------------------
TRIBITS_PACKAGE(salmon)
@@ -14,6 +17,9 @@ SET(HEADERS ExperimentalParameters.h
MultipleScattering.h
interface/cpp/MultipleScatteringInterface.h
+
+ Yield.h
+ interface/cpp/YieldInterface.h
)
APPEND_SET(SALMON_SOURCES
@@ -27,6 +33,11 @@ APPEND_SET(SALMON_SOURCES
interface/cpp/MultipleScatteringInterface.cpp
interface/fortran/MultipleScattering_M.f90
interface/fortran/MultipleScattering_I.f90
+
+ Yield.cpp
+ interface/cpp/YieldInterface.cpp
+ interface/fortran/Yield_M.f90
+ interface/fortran/Yield_I.f90
)
diff --git a/sammy/src/salmon/ExperimentalParameters.cpp b/sammy/src/salmon/ExperimentalParameters.cpp
index b1f6a370c519f5cb257d4a95db3fac08f05f9a51..5c0ed733e72fab5b73c48b69c00e98939a435ab1 100644
--- a/sammy/src/salmon/ExperimentalParameters.cpp
+++ b/sammy/src/salmon/ExperimentalParameters.cpp
@@ -8,6 +8,7 @@ namespace sammy{
radius = 0.0;
height = 0.0;
thickness = 0.0;
+ arealDensity = 0.0;
beamRadius = 0.0;
}
@@ -15,6 +16,7 @@ namespace sammy{
radius = orig.radius;
height = orig.height;
thickness = orig.thickness;
+ arealDensity = orig.arealDensity;
beamRadius = orig.beamRadius;
}
diff --git a/sammy/src/salmon/ExperimentalParameters.h b/sammy/src/salmon/ExperimentalParameters.h
index 59a60b5240dd5b6ac4df5bf1fbb9b7903c3f19d0..17eed0c4726a55d6590933ae4713d5b24225c7b1 100644
--- a/sammy/src/salmon/ExperimentalParameters.h
+++ b/sammy/src/salmon/ExperimentalParameters.h
@@ -35,6 +35,13 @@ namespace sammy{
void setThickness(double t){
thickness = t;
}
+ /**
+ * Set the areal density
+ * @param ad the areal density of the sample [at/b]
+ */
+ void setArealDensity(double ad){
+ arealDensity = ad;
+ }
/**
* Set the beam radius
* @param br the beam radius [cm]
@@ -72,6 +79,12 @@ namespace sammy{
double getThickness() const{
return thickness;
}
+ /**
+ * Get the thickness of the sample [at/b]
+ */
+ double getArealDensity() const{
+ return arealDensity;
+ }
/**
* Get the beam radius [cm]
*/
@@ -106,6 +119,11 @@ namespace sammy{
******************************/
double thickness; /**< thickness [cm] */
+ /******************************
+ * The thickness of the sample in atom/b
+ ******************************/
+ double arealDensity; /**< thickness [atom/b] */
+
/******************************
* The radius of the beam striking the sample in cm
******************************/
diff --git a/sammy/src/salmon/MultipleScattering.cpp b/sammy/src/salmon/MultipleScattering.cpp
index 0257a6664c92933b0d59d366f414de5c12d9aa27..c8332226437ce80d8fe388c1b2466417cb872039 100644
--- a/sammy/src/salmon/MultipleScattering.cpp
+++ b/sammy/src/salmon/MultipleScattering.cpp
@@ -5,6 +5,7 @@
namespace sammy{
MultipleScattering::MultipleScattering(){
+ numScatters = 0;
numTheta = 0;
numZHatIntegration = 0;
numBeamAreaIntegration = 0;
@@ -18,6 +19,7 @@ namespace sammy{
}
MultipleScattering::MultipleScattering(const MultipleScattering &orig){
+ numScatters = orig.numScatters;
numTheta = orig.numTheta;
numZHatIntegration = orig.numZHatIntegration;
numBeamAreaIntegration = orig.numBeamAreaIntegration;
diff --git a/sammy/src/salmon/MultipleScattering.h b/sammy/src/salmon/MultipleScattering.h
index 7d6d18ce5aeffb11329c4e7b4d432a20385b8b28..36fd5b95017011b096b02a9f91fe31f2a8fbb2af 100644
--- a/sammy/src/salmon/MultipleScattering.h
+++ b/sammy/src/salmon/MultipleScattering.h
@@ -14,54 +14,59 @@ namespace sammy{
* Setters
******************************/
+ /**
+ * Set number of scatters considered before neutron absorption
+ * @param numScat the number of scatters considered before neutron absorption
+ */
+ void setNumScatters(int numScat) { numScatters = numScat;}
/**
* Set number of integration points over theta
- * @ param Ntheta number of integration points over theta
+ * @param Ntheta number of integration points over theta
*/
void setNumTheta(int Ntheta) { numTheta = Ntheta;}
/**
* Set number of integration points in the z-direction
- * @ param Ngausz number of integration points in the z-direction (4, 8, or 16 Default = 16)
+ * @param Ngausz number of integration points in the z-direction (4, 8, or 16 Default = 16)
*/
void setNumZHatIntegration(int Ngausz) { numZHatIntegration = Ngausz;}
/**
* Set number of integration points over cross section of beam
- * @ param Ngaus number of integration points over cross section of beam (4, 8, or 16 Default = 16)
+ * @param Ngaus number of integration points over cross section of beam (4, 8, or 16 Default = 16)
*/
void setNumBeamAreaIntegration(int Ngaus) { numBeamAreaIntegration = Ngaus;}
/**
* Set number of theta integration points near the end points
- * @ param Mtheta number of integration points theta integration near the end points
+ * @param Mtheta number of integration points theta integration near the end points
*/
void setNumThetaEndPoints(int Mtheta) { numThetaEndPoints = Mtheta;}
/**
* Set number of integration points for interpolation of sigma
- * @ param Nxtpt number of integration points for interpolation of sigma
+ * @param Nxtpt number of integration points for interpolation of sigma
*/
void setNumInterpSigma(int Nxtpt) { numInterpSigma = Nxtpt;}
/**
* Set number of integration points for interpolation over sigma prime (scattered)
- * @ param Mxtpt number of integration points for interpolation over sigma prime (scattered)
+ * @param Mxtpt number of integration points for interpolation over sigma prime (scattered)
*/
void setNumInterpSigmaPrime(int Mxtpt) { numInterpSigmaPrime = Mxtpt;}
/**
* Set minimum value of log of total cross section
- * @ param Xtmin minimum value of log of total cross section
+ * @param Xtmin minimum value of log of total cross section
*/
void setLogSigmaTotMin(double Xtmin) { logSigmaTotMin = Xtmin;}
/**
* Set maximum value of log of total cross section
- * @ param Xtmax maximum value of log of total cross section
+ * @param Xtmax maximum value of log of total cross section
*/
void setLogSigmaTotMax(double Xtmax) { logSigmaTotMax = Xtmax;}
/**
* Set number of integration points in theta integration near cos(theta) = 1
- * @ param Jtheta number of integration points in theta integration near cos(theta) = 1
+ * @param Jtheta number of integration points in theta integration near cos(theta) = 1
*/
void setNumThetaNearOne(int Jtheta) { numThetaNearOne = Jtheta;}
/**
* Set number of integration points in theta integration near cos(theta) = 0
- * @ param Ktheta number of integration points in theta integration near cos(theta) = 0
+ * @param Ktheta number of integration points in theta integration near cos(theta) = 0
*/
void setNumThetaNearZero(int Ktheta) { numThetaNearZero = Ktheta;}
@@ -69,6 +74,12 @@ namespace sammy{
Getters
******************************/
+ /**
+ * Get number of scatters considered before neutron absorption
+ */
+ int getNumScatters() const{
+ return numScatters;
+ }
/**
* Get number of integration points over theta
*/
@@ -134,6 +145,11 @@ namespace sammy{
private:
+ /**
+ * How many scatters to consider before neutron absorption?
+ */
+ int numScatters;
+
/******************************
* The number of points to integrate over angle theta (Default = 33)
******************************/
diff --git a/sammy/src/salmon/Yield.cpp b/sammy/src/salmon/Yield.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..88e4588dbc5707063a222b035a6dbcf7004b8902
--- /dev/null
+++ b/sammy/src/salmon/Yield.cpp
@@ -0,0 +1,32 @@
+#include "Yield.h"
+#include <iostream>
+
+namespace sammy{
+
+ Yield::Yield()
+ {
+ normString = "";
+ normType = normalization::NOTSET;
+ calcDerivs = false;
+ numParam = 0;
+ Y0 = 0.0;
+ Y1 = 0.0;
+ Y2 = 0.0;
+ totY = 0.0;
+ }
+
+ Yield::Yield(const Yield &orig)
+ {
+ normType = orig.normType;
+ calcDerivs = orig.calcDerivs;
+ numParam = orig.numParam;
+ Y0 = orig.Y0;
+ Y1 = orig.Y1;
+ Y2 = orig.Y2;
+ totY = orig.totY;
+ derivY0.insert( derivY0.begin(), orig.derivY0.begin(), orig.derivY0.end() );
+ derivY1.insert( derivY1.begin(), orig.derivY1.begin(), orig.derivY1.end() );
+ derivY2.insert( derivY2.begin(), orig.derivY2.begin(), orig.derivY2.end() );
+ derivTotY.insert( derivTotY.begin(), orig.derivTotY.begin(), orig.derivTotY.end() );
+ }
+}
\ No newline at end of file
diff --git a/sammy/src/salmon/Yield.h b/sammy/src/salmon/Yield.h
new file mode 100644
index 0000000000000000000000000000000000000000..8eca77bfd0059c38471b865f0b2891ad17617e86
--- /dev/null
+++ b/sammy/src/salmon/Yield.h
@@ -0,0 +1,287 @@
+#ifndef YIELD_H
+#define YIELD_H
+
+#include <vector>
+#include <string>
+#include <memory>
+
+namespace sammy{
+
+ /** @brief Storage for the observed quantity of yield.
+
+ Store the values specific to yields (capture, fission, scattering...)
+ @date February 2020
+ */
+ class Yield{
+
+ public:
+
+
+ /******************************
+ * The type of normalization applied to the yield
+ ******************************/
+ enum normalization {SELFSHIELDED, DIVIDEBYN, TIMESSIGTOT, NOTSET} normType;
+
+ Yield();
+ Yield(const Yield& orig);
+ virtual ~Yield(){};
+
+ /******************************
+ * Setters
+ ******************************/
+
+ /**
+ * Set normalization type for yield to: self-shielded
+ */
+ void setNormTypeSelfShielded() {
+ normType = normalization::SELFSHIELDED;
+ normString = "Yield = Normalization type";
+ }
+
+ /**
+ * Set normalization type for yield to: divide by sample thickness \f$ N \f$
+ */
+ void setNormTypeDivideByN() {
+ normType = normalization::DIVIDEBYN;
+ normString = "Cross section = Normalization type";
+ }
+
+ /**
+ * Set normalization type for yield to: multiply by \f$ \sigma_t \f$
+ */
+ void setNormTypeTimesSigTot() {
+ normType = normalization::TIMESSIGTOT;
+ normString = "Other (1-e)sigma = Normalization type";
+ }
+
+ /**
+ * Set the boolean for whether to solve Bayes (calc. derivs.)
+ * @param sb boolean for whether we are solving Bayes eqtns.
+ */
+ void setCalcDerivs(bool cd) { calcDerivs = cd;}
+
+ /**
+ * Set the number of varied parameters (this should be temporary...we should take
+ * the info from a CovarianceData object or resonance parameter object)
+ */
+ void setNumTotalParam(int np){
+ numParam = np;
+ }
+
+ /**
+ * Set the primary yield
+ * @param y0 the primary yield
+ */
+ void setY0(double y0) { Y0 = y0;}
+
+ /**
+ * Set the secondary yield
+ * @param y1 the secondary yield
+ */
+ void setY1(double y1) { Y1 = y1;}
+
+ /**
+ * Set the tertiary yield
+ * @param y2 the tertiary yield
+ */
+ void setY2(double y2) { Y2 = y2;}
+
+ /**
+ * Set the total yield
+ * @param yt the total yield
+ */
+ void setTotY(double yt) { totY = yt; }
+
+ /**
+ * Set the derivative of primary yield w.r.t. parameter at indice iPar
+ * @param iPar index for the parameter of which deriv. was taken w.r.t.
+ * @param dy0 the primary yield
+ */
+ void setDerivY0(int iPar, double dy0){
+ if ((int)derivY0.size() < (iPar+1)){
+ derivY0.resize(iPar+1);
+ }
+ derivY0[iPar] = dy0;
+ }
+
+ /**
+ * Set the derivative of secondary yield w.r.t. parameter at indice iPar
+ * @param iPar index for the parameter of which deriv. was taken w.r.t.
+ * @param dy1 the secondary yield
+ */
+ void setDerivY1(int iPar, double dy1){
+ if ((int)derivY1.size() < (iPar+1)){
+ derivY1.resize(iPar+1);
+ }
+ derivY1[iPar] = dy1;
+ }
+
+ /**
+ * Set the derivative of tertiary yield w.r.t. parameter at indice iPar
+ * @param iPar index for the parameter of which deriv. was taken w.r.t.
+ * @param dy2 the tertiary yield
+ */
+ void setDerivY2(int iPar, double dy2){
+ if ((int)derivY2.size() < (iPar+1)){
+ derivY2.resize(iPar+1);
+ }
+ derivY2[iPar] = dy2;
+ }
+
+ /**
+ * Set the derivative of total yield w.r.t. parameter at indice iPar
+ * @param iPar index for the parameter of which deriv. was taken w.r.t.
+ * @param dyt the total yield
+ */
+ void setDerivTotY(int iPar, double dyt){
+ if ((int)derivTotY.size() < (iPar+1)){
+ derivTotY.resize(iPar+1);
+ }
+ derivTotY[iPar] = dyt;
+ }
+
+ /******************************
+ Getters
+ ******************************/
+
+ /**
+ * Get the normalization string for printing
+ */
+ std::string getNormString() const{ return normString; }
+
+ /**
+ * Get boolean whether normalization type for yield is self-shielded
+ */
+ bool normTypeIsSelfShielded() const{
+ return normType == normalization::SELFSHIELDED;
+ }
+
+ /**
+ * Get boolean whether normalization type for yield is divide by sample thickness \f$ N \f$
+ */
+ bool normTypeIsDivideByN() const{
+ return normType == normalization::DIVIDEBYN;
+ }
+
+ /**
+ * Get boolean whether normalization type for yield is to multiply by \f$ \sigma_t \f$
+ */
+ bool normTypeIsTimesSigTot() const{
+ return normType == normalization::TIMESSIGTOT;
+ }
+
+ /**
+ * Get boolean whether normalization type for yield has not been set
+ */
+ bool normTypeIsNotSet() const{
+ return normType == normalization::NOTSET;
+ }
+
+ /**
+ * Get decision on whether we solve Bayes (calc. derivs.)
+ */
+ bool getCalcDerivs() const{
+ return calcDerivs;
+ }
+
+ /**
+ * Get the total number of varied parameters (this should be temporary...we want to take
+ * the info from a CovarianceData object or resonance parameter object)
+ */
+ int getNumTotalParam() const{
+ return numParam;
+ }
+
+ /**
+ * Get the primary yield
+ */
+ double getY0() const{
+ return Y0;
+ }
+
+ /**
+ * Get the secondary yield
+ */
+ double getY1() const{
+ return Y1;
+ }
+
+ /**
+ * Get the tertiary yield
+ */
+ double getY2() const{
+ return Y2;
+ }
+
+ /**
+ * Get the total yield
+ */
+ double getTotY() const{
+ return totY;
+ }
+
+ /**
+ * Get the derivative of primary yield w.r.t. param. iPar
+ * @param iPar the indice of the parameter of which der. was taken w.r.t.
+ */
+ double getDerivY0(int iPar) const{
+ return derivY0[iPar];
+ }
+
+ /**
+ * Get the derivative of secondary yield w.r.t. param. iPar
+ * @param iPar the indice of the parameter of which der. was taken w.r.t.
+ */
+ double getDerivY1(int iPar) const{
+ return derivY1[iPar];
+ }
+
+ /**
+ * Get the derivative of tertiary yield w.r.t. param. iPar
+ * @param iPar the indice of the parameter of which der. was taken w.r.t.
+ */
+ double getDerivY2(int iPar) const{
+ return derivY2[iPar];
+ }
+
+ /**
+ * Get the derivative oftotal yield
+ */
+ double getDerivTotY(int iPar) const{
+ return derivTotY[iPar];
+ }
+
+ private:
+
+ /******************************
+ * String to print to files, stdout, etc.
+ *****************************/
+
+ std::string normString;
+
+ /******************************
+ * Should we calculate derivatives?
+ *****************************/
+
+ bool calcDerivs;
+
+ /******************************
+ * How many varied parameters are there?
+ *****************************/
+ int numParam;
+
+ /******************************
+ * The primary, secondary, tertiary, and total (sum) yields
+ ******************************/
+ double Y0, Y1, Y2, totY;
+
+ /******************************
+ * The derivatives of the primary, secondary, tertiary, and total (sum) yields
+ ******************************/
+ std::vector<double> derivY0, derivY1, derivY2, derivTotY;
+ };
+
+}
+
+
+#endif // YIELD_H
diff --git a/sammy/src/salmon/interface/cix/ExperimentalParameters.cpp2f.xml b/sammy/src/salmon/interface/cix/ExperimentalParameters.cpp2f.xml
index 4169973e69a9fe3861c6d5fa30db8a418d831418..40d773d09b34eecc00f8b755cf3cc5c2d29548c9 100644
--- a/sammy/src/salmon/interface/cix/ExperimentalParameters.cpp2f.xml
+++ b/sammy/src/salmon/interface/cix/ExperimentalParameters.cpp2f.xml
@@ -13,6 +13,9 @@
<method name="setThickness">
<param name="t" type="double"/>
</method>
+ <method name="setArealDensity">
+ <param name="t" type="double"/>
+ </method>
<method name="setBeamRadius">
<param name="br" type="double"/>
</method>
@@ -23,6 +26,7 @@
<method name="getRadius" return_type="double"/>
<method name="getHeight" return_type="double"/>
<method name="getThickness" return_type="double"/>
+ <method name="getArealDensity" return_type="double"/>
<method name="getBeamRadius" return_type="double"/>
<method name="getBeamHeight" return_type="double"/>
diff --git a/sammy/src/salmon/interface/cix/MultipleScattering.cpp2f.xml b/sammy/src/salmon/interface/cix/MultipleScattering.cpp2f.xml
index ae7d1a39854184fe221c712fada24fa3be03261f..6262f21a1be3cfe21968fe9514a9a8e4ebdd24e5 100644
--- a/sammy/src/salmon/interface/cix/MultipleScattering.cpp2f.xml
+++ b/sammy/src/salmon/interface/cix/MultipleScattering.cpp2f.xml
@@ -1,9 +1,12 @@
<generate name="MultipleScattering">
- <include_relative name="../MultipleScattering.h"/>
+ <include_relative name="../../MultipleScattering.h"/>
<using_namespace name="sammy"/>
<class name="MultipleScattering">
+ <method name="setNumScatters">
+ <param name="ns" type="int"/>
+ </method>
<method name="setNumTheta">
<param name="nt" type="int"/>
</method>
@@ -35,6 +38,7 @@
<param name="ntnearzero" type="int"/>
</method>
+ <method name="getNumScatters" return_type="int"/>
<method name="getNumTheta" return_type="int"/>
<method name="getNumZHatIntegration" return_type="int"/>
<method name="getNumBeamAreaIntegration" return_type="int"/>
@@ -47,4 +51,4 @@
<method name="getNumThetaNearZero" return_type="int"/>
</class>
-</generate>
\ No newline at end of file
+</generate>
diff --git a/sammy/src/salmon/interface/cix/Yield.cpp2f.xml b/sammy/src/salmon/interface/cix/Yield.cpp2f.xml
new file mode 100644
index 0000000000000000000000000000000000000000..1fe8432955162291eacf8a3b28ca09cf88bdffc6
--- /dev/null
+++ b/sammy/src/salmon/interface/cix/Yield.cpp2f.xml
@@ -0,0 +1,70 @@
+<generate name="Yield">
+ <include_relative name="../../Yield.h"/>
+ <using_namespace name="sammy"/>
+
+ <class name="Yield">
+
+ <method name="setNormTypeSelfShielded"/>
+ <method name="setNormTypeDivideByN"/>
+ <method name="setNormTypeTimesSigTot"/>
+ <method name="setCalcDerivs">
+ <param name="sb" type="bool"/>
+ </method>
+ <method name="setNumTotalParam">
+ <param name="np" type="int"/>
+ </method>
+ <method name="setY0">
+ <param name="y0" type="double"/>
+ </method>
+ <method name="setY1">
+ <param name="y1" type="double"/>
+ </method>
+ <method name="setY2">
+ <param name="y2" type="double"/>
+ </method>
+ <method name="setTotY">
+ <param name="toty" type="double"/>
+ </method>
+ <method name="setDerivY0">
+ <param name="iPar" type="int" offset="-1"/>
+ <param name="dy0" type="double"/>
+ </method>
+ <method name="setDerivY1">
+ <param name="iPar" type="int" offset="-1"/>
+ <param name="dy1" type="double"/>
+ </method>
+ <method name="setDerivY2">
+ <param name="iPar" type="int" offset="-1"/>
+ <param name="dy2" type="double"/>
+ </method>
+ <method name="setDerivTotY">
+ <param name="iPar" type="int" offset="-1"/>
+ <param name="dy2" type="double"/>
+ </method>
+
+ <method name="getNormString" return_type="string"/>
+ <method name="normTypeIsSelfShielded" return_type="bool"/>
+ <method name="normTypeIsDivideByN" return_type="bool"/>
+ <method name="normTypeIsTimesSigTot" return_type="bool"/>
+ <method name="getCalcDerivs" return_type="bool"/>
+ <method name="getNumTotalParam" return_type="int"/>
+ <method name="getY0" return_type="double"/>
+ <method name="getY1" return_type="double"/>
+ <method name="getY2" return_type="double"/>
+ <method name="getTotY" return_type="double"/>
+
+ <method name="getDerivY0" return_type="double">
+ <param name="iPar" type="int" offset="-1"/>
+ </method>
+ <method name="getDerivY1" return_type="double">
+ <param name="iPar" type="int" offset="-1"/>
+ </method>
+ <method name="getDerivY2" return_type="double">
+ <param name="iPar" type="int" offset="-1"/>
+ </method>
+ <method name="getDerivTotY" return_type="double">
+ <param name="iPar" type="int" offset="-1"/>
+ </method>
+
+ </class>
+</generate>
\ No newline at end of file
diff --git a/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.cpp b/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.cpp
index c3eacbacf866d2a17538326c1d20ba5a075b7c51..96e23ad4c73b647548c9c1f2900ddae5d466d4db 100644
--- a/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.cpp
+++ b/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.cpp
@@ -2,7 +2,7 @@
* This file has been dynamically generated by Class Interface Xml (CIX)
* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
* If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Fri Dec 20 13:08:25 EST 2019
+* Date Generated: Tue Mar 03 09:03:32 EST 2020
* If any issues are experiences with this generated file that cannot be fixed
* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
*/
@@ -24,6 +24,11 @@ void SampleDimensions_setThickness(void * SampleDimensions_ptr,double * t)
((SampleDimensions*)SampleDimensions_ptr)->setThickness(*t);
}
+void SampleDimensions_setArealDensity(void * SampleDimensions_ptr,double * t)
+{
+ ((SampleDimensions*)SampleDimensions_ptr)->setArealDensity(*t);
+}
+
void SampleDimensions_setBeamRadius(void * SampleDimensions_ptr,double * br)
{
((SampleDimensions*)SampleDimensions_ptr)->setBeamRadius(*br);
@@ -49,6 +54,11 @@ double SampleDimensions_getThickness(void * SampleDimensions_ptr)
return ((SampleDimensions*)SampleDimensions_ptr)->getThickness();
}
+double SampleDimensions_getArealDensity(void * SampleDimensions_ptr)
+{
+ return ((SampleDimensions*)SampleDimensions_ptr)->getArealDensity();
+}
+
double SampleDimensions_getBeamRadius(void * SampleDimensions_ptr)
{
return ((SampleDimensions*)SampleDimensions_ptr)->getBeamRadius();
diff --git a/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.h b/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.h
index a33085c0eb4a4dce1e06eb70cda596bfe544476d..43d5407dc7d1bdfe7d8241db9019d702b3d0f647 100644
--- a/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.h
+++ b/sammy/src/salmon/interface/cpp/ExperimentalParametersInterface.h
@@ -2,7 +2,7 @@
* This file has been dynamically generated by Class Interface Xml (CIX)
* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
* If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Fri Dec 20 13:08:25 EST 2019
+* Date Generated: Tue Mar 03 09:03:32 EST 2020
* If any issues are experiences with this generated file that cannot be fixed
* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
*/
@@ -16,11 +16,13 @@ extern "C" {
void SampleDimensions_setRadius(void * SampleDimensions_ptr,double * r);
void SampleDimensions_setHeight(void * SampleDimensions_ptr,double * h);
void SampleDimensions_setThickness(void * SampleDimensions_ptr,double * t);
+void SampleDimensions_setArealDensity(void * SampleDimensions_ptr,double * t);
void SampleDimensions_setBeamRadius(void * SampleDimensions_ptr,double * br);
void SampleDimensions_setBeamHeight(void * SampleDimensions_ptr,double * br);
double SampleDimensions_getRadius(void * SampleDimensions_ptr);
double SampleDimensions_getHeight(void * SampleDimensions_ptr);
double SampleDimensions_getThickness(void * SampleDimensions_ptr);
+double SampleDimensions_getArealDensity(void * SampleDimensions_ptr);
double SampleDimensions_getBeamRadius(void * SampleDimensions_ptr);
double SampleDimensions_getBeamHeight(void * SampleDimensions_ptr);
void* SampleDimensions_initialize();
diff --git a/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.cpp b/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.cpp
index 1e88a180d7713e116d08312a00a09514588ed95b..65c95b66f5590248f1edea13118ead03fe1154ce 100644
--- a/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.cpp
+++ b/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.cpp
@@ -2,13 +2,18 @@
* This file has been dynamically generated by Class Interface Xml (CIX)
* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
* If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Mon Dec 09 09:48:43 EST 2019
+* Date Generated: Mon Mar 02 10:33:41 EST 2020
* If any issues are experiences with this generated file that cannot be fixed
* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
*/
#include <string.h>
#include "MultipleScatteringInterface.h"
using namespace sammy;
+void MultipleScattering_setNumScatters(void * MultipleScattering_ptr,int * ns)
+{
+ ((MultipleScattering*)MultipleScattering_ptr)->setNumScatters(*ns);
+}
+
void MultipleScattering_setNumTheta(void * MultipleScattering_ptr,int * nt)
{
((MultipleScattering*)MultipleScattering_ptr)->setNumTheta(*nt);
@@ -59,6 +64,11 @@ void MultipleScattering_setNumThetaNearZero(void * MultipleScattering_ptr,int *
((MultipleScattering*)MultipleScattering_ptr)->setNumThetaNearZero(*ntnearzero);
}
+int MultipleScattering_getNumScatters(void * MultipleScattering_ptr)
+{
+ return ((MultipleScattering*)MultipleScattering_ptr)->getNumScatters();
+}
+
int MultipleScattering_getNumTheta(void * MultipleScattering_ptr)
{
return ((MultipleScattering*)MultipleScattering_ptr)->getNumTheta();
diff --git a/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.h b/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.h
index 269bd916d3f88b268bb0b69a949ed92432618a7f..778c51f0677a1a58281eccd04759c3cd2832304c 100644
--- a/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.h
+++ b/sammy/src/salmon/interface/cpp/MultipleScatteringInterface.h
@@ -2,7 +2,7 @@
* This file has been dynamically generated by Class Interface Xml (CIX)
* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
* If changes need to occur, modify the appropriate CIX xml file
-* Date Generated: Mon Dec 09 09:48:43 EST 2019
+* Date Generated: Mon Mar 02 10:33:41 EST 2020
* If any issues are experiences with this generated file that cannot be fixed
* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
*/
@@ -13,6 +13,7 @@ using namespace sammy;
#ifdef __cplusplus
extern "C" {
#endif
+void MultipleScattering_setNumScatters(void * MultipleScattering_ptr,int * ns);
void MultipleScattering_setNumTheta(void * MultipleScattering_ptr,int * nt);
void MultipleScattering_setNumZHatIntegration(void * MultipleScattering_ptr,int * nz);
void MultipleScattering_setNumBeamAreaIntegration(void * MultipleScattering_ptr,int * nba);
@@ -23,6 +24,7 @@ void MultipleScattering_setLogSigmaTotMin(void * MultipleScattering_ptr,double *
void MultipleScattering_setLogSigmaTotMax(void * MultipleScattering_ptr,double * nsigmax);
void MultipleScattering_setNumThetaNearOne(void * MultipleScattering_ptr,int * ntnearone);
void MultipleScattering_setNumThetaNearZero(void * MultipleScattering_ptr,int * ntnearzero);
+int MultipleScattering_getNumScatters(void * MultipleScattering_ptr);
int MultipleScattering_getNumTheta(void * MultipleScattering_ptr);
int MultipleScattering_getNumZHatIntegration(void * MultipleScattering_ptr);
int MultipleScattering_getNumBeamAreaIntegration(void * MultipleScattering_ptr);
diff --git a/sammy/src/salmon/interface/cpp/YieldInterface.cpp b/sammy/src/salmon/interface/cpp/YieldInterface.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..e118efb1b69aeb0281e6e74b448b40ba05797167
--- /dev/null
+++ b/sammy/src/salmon/interface/cpp/YieldInterface.cpp
@@ -0,0 +1,163 @@
+/*!
+* This file has been dynamically generated by Class Interface Xml (CIX)
+* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+* If changes need to occur, modify the appropriate CIX xml file
+* Date Generated: Wed Apr 08 13:29:36 EDT 2020
+* If any issues are experiences with this generated file that cannot be fixed
+* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+*/
+#include <string.h>
+#include "YieldInterface.h"
+using namespace sammy;
+void Yield_setNormTypeSelfShielded(void * Yield_ptr)
+{
+ ((Yield*)Yield_ptr)->setNormTypeSelfShielded();
+}
+
+void Yield_setNormTypeDivideByN(void * Yield_ptr)
+{
+ ((Yield*)Yield_ptr)->setNormTypeDivideByN();
+}
+
+void Yield_setNormTypeTimesSigTot(void * Yield_ptr)
+{
+ ((Yield*)Yield_ptr)->setNormTypeTimesSigTot();
+}
+
+void Yield_setCalcDerivs(void * Yield_ptr,bool * sb)
+{
+ ((Yield*)Yield_ptr)->setCalcDerivs(*sb);
+}
+
+void Yield_setNumTotalParam(void * Yield_ptr,int * np)
+{
+ ((Yield*)Yield_ptr)->setNumTotalParam(*np);
+}
+
+void Yield_setY0(void * Yield_ptr,double * y0)
+{
+ ((Yield*)Yield_ptr)->setY0(*y0);
+}
+
+void Yield_setY1(void * Yield_ptr,double * y1)
+{
+ ((Yield*)Yield_ptr)->setY1(*y1);
+}
+
+void Yield_setY2(void * Yield_ptr,double * y2)
+{
+ ((Yield*)Yield_ptr)->setY2(*y2);
+}
+
+void Yield_setTotY(void * Yield_ptr,double * toty)
+{
+ ((Yield*)Yield_ptr)->setTotY(*toty);
+}
+
+void Yield_setDerivY0(void * Yield_ptr,int * iPar,double * dy0)
+{
+ ((Yield*)Yield_ptr)->setDerivY0(*iPar,*dy0);
+}
+
+void Yield_setDerivY1(void * Yield_ptr,int * iPar,double * dy1)
+{
+ ((Yield*)Yield_ptr)->setDerivY1(*iPar,*dy1);
+}
+
+void Yield_setDerivY2(void * Yield_ptr,int * iPar,double * dy2)
+{
+ ((Yield*)Yield_ptr)->setDerivY2(*iPar,*dy2);
+}
+
+void Yield_setDerivTotY(void * Yield_ptr,int * iPar,double * dy2)
+{
+ ((Yield*)Yield_ptr)->setDerivTotY(*iPar,*dy2);
+}
+
+void Yield_getNormString(void * Yield_ptr, char * string2Populate, int * lengthOfString2Populate)
+{
+ const std::string& cix_string = ((Yield*)Yield_ptr)->getNormString();
+ int t_length = *lengthOfString2Populate < static_cast<int>(cix_string.length())
+ ? *lengthOfString2Populate : cix_string.length();
+ strncpy(string2Populate, cix_string.c_str(), t_length);
+ for(int cix_i = t_length; cix_i < *lengthOfString2Populate; cix_i++){
+ string2Populate[cix_i] = ' ';
+ }
+
+}
+
+bool Yield_normTypeIsSelfShielded(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->normTypeIsSelfShielded();
+}
+
+bool Yield_normTypeIsDivideByN(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->normTypeIsDivideByN();
+}
+
+bool Yield_normTypeIsTimesSigTot(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->normTypeIsTimesSigTot();
+}
+
+bool Yield_getCalcDerivs(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getCalcDerivs();
+}
+
+int Yield_getNumTotalParam(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getNumTotalParam();
+}
+
+double Yield_getY0(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getY0();
+}
+
+double Yield_getY1(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getY1();
+}
+
+double Yield_getY2(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getY2();
+}
+
+double Yield_getTotY(void * Yield_ptr)
+{
+ return ((Yield*)Yield_ptr)->getTotY();
+}
+
+double Yield_getDerivY0(void * Yield_ptr,int * iPar)
+{
+ return ((Yield*)Yield_ptr)->getDerivY0(*iPar);
+}
+
+double Yield_getDerivY1(void * Yield_ptr,int * iPar)
+{
+ return ((Yield*)Yield_ptr)->getDerivY1(*iPar);
+}
+
+double Yield_getDerivY2(void * Yield_ptr,int * iPar)
+{
+ return ((Yield*)Yield_ptr)->getDerivY2(*iPar);
+}
+
+double Yield_getDerivTotY(void * Yield_ptr,int * iPar)
+{
+ return ((Yield*)Yield_ptr)->getDerivTotY(*iPar);
+}
+
+void* Yield_initialize()
+{
+ return new Yield();
+}
+
+void Yield_destroy(void * Yield_ptr)
+{
+ delete (Yield*)Yield_ptr;
+}
+
diff --git a/sammy/src/salmon/interface/cpp/YieldInterface.h b/sammy/src/salmon/interface/cpp/YieldInterface.h
new file mode 100644
index 0000000000000000000000000000000000000000..64065a15ad56d7345537b199f7cdb1012d706a6e
--- /dev/null
+++ b/sammy/src/salmon/interface/cpp/YieldInterface.h
@@ -0,0 +1,48 @@
+/*!
+* This file has been dynamically generated by Class Interface Xml (CIX)
+* DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+* If changes need to occur, modify the appropriate CIX xml file
+* Date Generated: Wed Apr 08 13:29:36 EDT 2020
+* If any issues are experiences with this generated file that cannot be fixed
+* with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+*/
+#ifndef YIELDINTERFACE_H
+#define YIELDINTERFACE_H
+#include "../../Yield.h"
+using namespace sammy;
+#ifdef __cplusplus
+extern "C" {
+#endif
+void Yield_setNormTypeSelfShielded(void * Yield_ptr);
+void Yield_setNormTypeDivideByN(void * Yield_ptr);
+void Yield_setNormTypeTimesSigTot(void * Yield_ptr);
+void Yield_setCalcDerivs(void * Yield_ptr,bool * sb);
+void Yield_setNumTotalParam(void * Yield_ptr,int * np);
+void Yield_setY0(void * Yield_ptr,double * y0);
+void Yield_setY1(void * Yield_ptr,double * y1);
+void Yield_setY2(void * Yield_ptr,double * y2);
+void Yield_setTotY(void * Yield_ptr,double * toty);
+void Yield_setDerivY0(void * Yield_ptr,int * iPar,double * dy0);
+void Yield_setDerivY1(void * Yield_ptr,int * iPar,double * dy1);
+void Yield_setDerivY2(void * Yield_ptr,int * iPar,double * dy2);
+void Yield_setDerivTotY(void * Yield_ptr,int * iPar,double * dy2);
+void Yield_getNormString(void * Yield_ptr, char * string2Populate, int * lengthOfString2Populate);
+bool Yield_normTypeIsSelfShielded(void * Yield_ptr);
+bool Yield_normTypeIsDivideByN(void * Yield_ptr);
+bool Yield_normTypeIsTimesSigTot(void * Yield_ptr);
+bool Yield_getCalcDerivs(void * Yield_ptr);
+int Yield_getNumTotalParam(void * Yield_ptr);
+double Yield_getY0(void * Yield_ptr);
+double Yield_getY1(void * Yield_ptr);
+double Yield_getY2(void * Yield_ptr);
+double Yield_getTotY(void * Yield_ptr);
+double Yield_getDerivY0(void * Yield_ptr,int * iPar);
+double Yield_getDerivY1(void * Yield_ptr,int * iPar);
+double Yield_getDerivY2(void * Yield_ptr,int * iPar);
+double Yield_getDerivTotY(void * Yield_ptr,int * iPar);
+void* Yield_initialize();
+void Yield_destroy(void * Yield_ptr);
+#ifdef __cplusplus
+}
+#endif
+#endif /* YIELDINTERFACE_H */
diff --git a/sammy/src/salmon/interface/fortran/ExperimentalParameters_I.f90 b/sammy/src/salmon/interface/fortran/ExperimentalParameters_I.f90
index 060172f1067e28727fc42df46822ae721c1d524a..f51a5731538530273344f3fea49301687ba65eb3 100644
--- a/sammy/src/salmon/interface/fortran/ExperimentalParameters_I.f90
+++ b/sammy/src/salmon/interface/fortran/ExperimentalParameters_I.f90
@@ -2,7 +2,7 @@
!! This file has been dynamically generated by Class Interface Xml (CIX)
!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
!! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Fri Dec 20 13:08:25 EST 2019
+!! Date Generated: Tue Mar 03 09:03:32 EST 2020
!! If any issues are experiences with this generated file that cannot be fixed
!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
!!/
@@ -27,6 +27,12 @@ subroutine f_SampleDimensions_setThickness(SampleDimensions_ptr, t ) BIND(C,name
type(C_PTR), value :: SampleDimensions_ptr;
real(C_DOUBLE) :: t;
end subroutine
+subroutine f_SampleDimensions_setArealDensity(SampleDimensions_ptr, t ) BIND(C,name="SampleDimensions_setArealDensity")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: SampleDimensions_ptr;
+ real(C_DOUBLE) :: t;
+end subroutine
subroutine f_SampleDimensions_setBeamRadius(SampleDimensions_ptr, br ) BIND(C,name="SampleDimensions_setBeamRadius")
use,intrinsic :: ISO_C_BINDING
implicit none
@@ -54,6 +60,11 @@ real(C_DOUBLE) function f_SampleDimensions_getThickness(SampleDimensions_ptr ) B
implicit none
type(C_PTR), value :: SampleDimensions_ptr;
end function
+real(C_DOUBLE) function f_SampleDimensions_getArealDensity(SampleDimensions_ptr ) BIND(C,name="SampleDimensions_getArealDensity")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: SampleDimensions_ptr;
+end function
real(C_DOUBLE) function f_SampleDimensions_getBeamRadius(SampleDimensions_ptr ) BIND(C,name="SampleDimensions_getBeamRadius")
use,intrinsic :: ISO_C_BINDING
implicit none
diff --git a/sammy/src/salmon/interface/fortran/ExperimentalParameters_M.f90 b/sammy/src/salmon/interface/fortran/ExperimentalParameters_M.f90
index b24a121f30cc8bfe49c84969dacee0614f006b7d..b608ed5de3018991757ebeb31cad5414e7f56f18 100644
--- a/sammy/src/salmon/interface/fortran/ExperimentalParameters_M.f90
+++ b/sammy/src/salmon/interface/fortran/ExperimentalParameters_M.f90
@@ -2,7 +2,7 @@
!! This file has been dynamically generated by Class Interface Xml (CIX)
!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
!! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Fri Dec 20 13:08:25 EST 2019
+!! Date Generated: Tue Mar 03 09:03:32 EST 2020
!! If any issues are experiences with this generated file that cannot be fixed
!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
!!/
@@ -15,11 +15,13 @@ type SampleDimensions
procedure, pass(this) :: setRadius => SampleDimensions_setRadius
procedure, pass(this) :: setHeight => SampleDimensions_setHeight
procedure, pass(this) :: setThickness => SampleDimensions_setThickness
+ procedure, pass(this) :: setArealDensity => SampleDimensions_setArealDensity
procedure, pass(this) :: setBeamRadius => SampleDimensions_setBeamRadius
procedure, pass(this) :: setBeamHeight => SampleDimensions_setBeamHeight
procedure, pass(this) :: getRadius => SampleDimensions_getRadius
procedure, pass(this) :: getHeight => SampleDimensions_getHeight
procedure, pass(this) :: getThickness => SampleDimensions_getThickness
+ procedure, pass(this) :: getArealDensity => SampleDimensions_getArealDensity
procedure, pass(this) :: getBeamRadius => SampleDimensions_getBeamRadius
procedure, pass(this) :: getBeamHeight => SampleDimensions_getBeamHeight
procedure, pass(this) :: initialize => SampleDimensions_initialize
@@ -54,6 +56,12 @@ subroutine SampleDimensions_setThickness(this, t)
real(C_DOUBLE)::t
call f_SampleDimensions_setThickness(this%instance_ptr, t)
end subroutine
+subroutine SampleDimensions_setArealDensity(this, t)
+ implicit none
+ class(SampleDimensions)::this
+ real(C_DOUBLE)::t
+ call f_SampleDimensions_setArealDensity(this%instance_ptr, t)
+end subroutine
subroutine SampleDimensions_setBeamRadius(this, br)
implicit none
class(SampleDimensions)::this
@@ -84,6 +92,12 @@ function SampleDimensions_getThickness(this) result(result2Return)
real(C_DOUBLE):: result2Return
result2Return=f_SampleDimensions_getThickness(this%instance_ptr)
end function
+function SampleDimensions_getArealDensity(this) result(result2Return)
+ implicit none
+ class(SampleDimensions)::this
+ real(C_DOUBLE):: result2Return
+ result2Return=f_SampleDimensions_getArealDensity(this%instance_ptr)
+end function
function SampleDimensions_getBeamRadius(this) result(result2Return)
implicit none
class(SampleDimensions)::this
diff --git a/sammy/src/salmon/interface/fortran/MultipleScattering_I.f90 b/sammy/src/salmon/interface/fortran/MultipleScattering_I.f90
index f25f04dc2b093a98a4c3a1ca2bb3f95db5079420..5422e5029b773fad5bd1e34ad3825749d1a76600 100644
--- a/sammy/src/salmon/interface/fortran/MultipleScattering_I.f90
+++ b/sammy/src/salmon/interface/fortran/MultipleScattering_I.f90
@@ -2,13 +2,19 @@
!! This file has been dynamically generated by Class Interface Xml (CIX)
!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
!! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Mon Dec 09 09:48:43 EST 2019
+!! Date Generated: Mon Mar 02 10:33:41 EST 2020
!! If any issues are experiences with this generated file that cannot be fixed
!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
!!/
module MultipleScattering_I
use, intrinsic :: ISO_C_BINDING
interface
+subroutine f_MultipleScattering_setNumScatters(MultipleScattering_ptr, ns ) BIND(C,name="MultipleScattering_setNumScatters")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: MultipleScattering_ptr;
+ integer(C_INT) :: ns;
+end subroutine
subroutine f_MultipleScattering_setNumTheta(MultipleScattering_ptr, nt ) BIND(C,name="MultipleScattering_setNumTheta")
use,intrinsic :: ISO_C_BINDING
implicit none
@@ -69,6 +75,11 @@ subroutine f_MultipleScattering_setNumThetaNearZero(MultipleScattering_ptr, ntne
type(C_PTR), value :: MultipleScattering_ptr;
integer(C_INT) :: ntnearzero;
end subroutine
+integer(C_INT) function f_MultipleScattering_getNumScatters(MultipleScattering_ptr ) BIND(C,name="MultipleScattering_getNumScatters")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: MultipleScattering_ptr;
+end function
integer(C_INT) function f_MultipleScattering_getNumTheta(MultipleScattering_ptr ) BIND(C,name="MultipleScattering_getNumTheta")
use,intrinsic :: ISO_C_BINDING
implicit none
diff --git a/sammy/src/salmon/interface/fortran/MultipleScattering_M.f90 b/sammy/src/salmon/interface/fortran/MultipleScattering_M.f90
index aa10aee1db1b3d341a91a9300853d87958056e12..d6ed8dd22ff38a38a35860d32bd8f668ab09b0f3 100644
--- a/sammy/src/salmon/interface/fortran/MultipleScattering_M.f90
+++ b/sammy/src/salmon/interface/fortran/MultipleScattering_M.f90
@@ -2,7 +2,7 @@
!! This file has been dynamically generated by Class Interface Xml (CIX)
!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
!! If changes need to occur, modify the appropriate CIX xml file
-!! Date Generated: Mon Dec 09 09:48:43 EST 2019
+!! Date Generated: Mon Mar 02 10:33:41 EST 2020
!! If any issues are experiences with this generated file that cannot be fixed
!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
!!/
@@ -12,6 +12,7 @@ use MultipleScattering_I
type MultipleScattering
type(C_PTR) :: instance_ptr=C_NULL_PTR
contains
+ procedure, pass(this) :: setNumScatters => MultipleScattering_setNumScatters
procedure, pass(this) :: setNumTheta => MultipleScattering_setNumTheta
procedure, pass(this) :: setNumZHatIntegration => MultipleScattering_setNumZHatIntegration
procedure, pass(this) :: setNumBeamAreaIntegration => MultipleScattering_setNumBeamAreaIntegration
@@ -22,6 +23,7 @@ type MultipleScattering
procedure, pass(this) :: setLogSigmaTotMax => MultipleScattering_setLogSigmaTotMax
procedure, pass(this) :: setNumThetaNearOne => MultipleScattering_setNumThetaNearOne
procedure, pass(this) :: setNumThetaNearZero => MultipleScattering_setNumThetaNearZero
+ procedure, pass(this) :: getNumScatters => MultipleScattering_getNumScatters
procedure, pass(this) :: getNumTheta => MultipleScattering_getNumTheta
procedure, pass(this) :: getNumZHatIntegration => MultipleScattering_getNumZHatIntegration
procedure, pass(this) :: getNumBeamAreaIntegration => MultipleScattering_getNumBeamAreaIntegration
@@ -36,6 +38,12 @@ type MultipleScattering
procedure, pass(this) :: destroy => MultipleScattering_destroy
end type MultipleScattering
contains
+subroutine MultipleScattering_setNumScatters(this, ns)
+ implicit none
+ class(MultipleScattering)::this
+ integer(C_INT)::ns
+ call f_MultipleScattering_setNumScatters(this%instance_ptr, ns)
+end subroutine
subroutine MultipleScattering_setNumTheta(this, nt)
implicit none
class(MultipleScattering)::this
@@ -96,6 +104,12 @@ subroutine MultipleScattering_setNumThetaNearZero(this, ntnearzero)
integer(C_INT)::ntnearzero
call f_MultipleScattering_setNumThetaNearZero(this%instance_ptr, ntnearzero)
end subroutine
+function MultipleScattering_getNumScatters(this) result(result2Return)
+ implicit none
+ class(MultipleScattering)::this
+ integer(C_INT):: result2Return
+ result2Return=f_MultipleScattering_getNumScatters(this%instance_ptr)
+end function
function MultipleScattering_getNumTheta(this) result(result2Return)
implicit none
class(MultipleScattering)::this
diff --git a/sammy/src/salmon/interface/fortran/Yield_I.f90 b/sammy/src/salmon/interface/fortran/Yield_I.f90
new file mode 100644
index 0000000000000000000000000000000000000000..695d2115f39d94b4acb50b9147e8651f31add04d
--- /dev/null
+++ b/sammy/src/salmon/interface/fortran/Yield_I.f90
@@ -0,0 +1,178 @@
+!>
+!! This file has been dynamically generated by Class Interface Xml (CIX)
+!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+!! If changes need to occur, modify the appropriate CIX xml file
+!! Date Generated: Wed Apr 08 13:29:36 EDT 2020
+!! If any issues are experiences with this generated file that cannot be fixed
+!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+!!/
+module Yield_I
+use, intrinsic :: ISO_C_BINDING
+interface
+subroutine f_Yield_setNormTypeSelfShielded(Yield_ptr ) BIND(C,name="Yield_setNormTypeSelfShielded")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end subroutine
+subroutine f_Yield_setNormTypeDivideByN(Yield_ptr ) BIND(C,name="Yield_setNormTypeDivideByN")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end subroutine
+subroutine f_Yield_setNormTypeTimesSigTot(Yield_ptr ) BIND(C,name="Yield_setNormTypeTimesSigTot")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end subroutine
+subroutine f_Yield_setCalcDerivs(Yield_ptr, sb ) BIND(C,name="Yield_setCalcDerivs")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ logical(C_BOOL) :: sb;
+end subroutine
+subroutine f_Yield_setNumTotalParam(Yield_ptr, np ) BIND(C,name="Yield_setNumTotalParam")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: np;
+end subroutine
+subroutine f_Yield_setY0(Yield_ptr, y0 ) BIND(C,name="Yield_setY0")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ real(C_DOUBLE) :: y0;
+end subroutine
+subroutine f_Yield_setY1(Yield_ptr, y1 ) BIND(C,name="Yield_setY1")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ real(C_DOUBLE) :: y1;
+end subroutine
+subroutine f_Yield_setY2(Yield_ptr, y2 ) BIND(C,name="Yield_setY2")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ real(C_DOUBLE) :: y2;
+end subroutine
+subroutine f_Yield_setTotY(Yield_ptr, toty ) BIND(C,name="Yield_setTotY")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ real(C_DOUBLE) :: toty;
+end subroutine
+subroutine f_Yield_setDerivY0(Yield_ptr, iPar,dy0 ) BIND(C,name="Yield_setDerivY0")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+ real(C_DOUBLE) :: dy0;
+end subroutine
+subroutine f_Yield_setDerivY1(Yield_ptr, iPar,dy1 ) BIND(C,name="Yield_setDerivY1")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+ real(C_DOUBLE) :: dy1;
+end subroutine
+subroutine f_Yield_setDerivY2(Yield_ptr, iPar,dy2 ) BIND(C,name="Yield_setDerivY2")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+ real(C_DOUBLE) :: dy2;
+end subroutine
+subroutine f_Yield_setDerivTotY(Yield_ptr, iPar,dy2 ) BIND(C,name="Yield_setDerivTotY")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+ real(C_DOUBLE) :: dy2;
+end subroutine
+subroutine f_Yield_getNormString(Yield_ptr, string2Populate, lengthOfString2Populate) BIND(C,name="Yield_getNormString")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ character(C_CHAR) :: string2Populate;
+ integer(C_INT) :: lengthOfString2Populate;
+end subroutine
+logical(C_BOOL) function f_Yield_normTypeIsSelfShielded(Yield_ptr ) BIND(C,name="Yield_normTypeIsSelfShielded")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+logical(C_BOOL) function f_Yield_normTypeIsDivideByN(Yield_ptr ) BIND(C,name="Yield_normTypeIsDivideByN")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+logical(C_BOOL) function f_Yield_normTypeIsTimesSigTot(Yield_ptr ) BIND(C,name="Yield_normTypeIsTimesSigTot")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+logical(C_BOOL) function f_Yield_getCalcDerivs(Yield_ptr ) BIND(C,name="Yield_getCalcDerivs")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+integer(C_INT) function f_Yield_getNumTotalParam(Yield_ptr ) BIND(C,name="Yield_getNumTotalParam")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+real(C_DOUBLE) function f_Yield_getY0(Yield_ptr ) BIND(C,name="Yield_getY0")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+real(C_DOUBLE) function f_Yield_getY1(Yield_ptr ) BIND(C,name="Yield_getY1")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+real(C_DOUBLE) function f_Yield_getY2(Yield_ptr ) BIND(C,name="Yield_getY2")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+real(C_DOUBLE) function f_Yield_getTotY(Yield_ptr ) BIND(C,name="Yield_getTotY")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+end function
+real(C_DOUBLE) function f_Yield_getDerivY0(Yield_ptr, iPar ) BIND(C,name="Yield_getDerivY0")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+end function
+real(C_DOUBLE) function f_Yield_getDerivY1(Yield_ptr, iPar ) BIND(C,name="Yield_getDerivY1")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+end function
+real(C_DOUBLE) function f_Yield_getDerivY2(Yield_ptr, iPar ) BIND(C,name="Yield_getDerivY2")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+end function
+real(C_DOUBLE) function f_Yield_getDerivTotY(Yield_ptr, iPar ) BIND(C,name="Yield_getDerivTotY")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: Yield_ptr;
+ integer(C_INT) :: iPar;
+end function
+type(C_PTR) function f_Yield_initialize( )BIND(C,name="Yield_initialize")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR) :: Yield_ptr;
+end function
+subroutine f_Yield_destroy(this) BIND(C,name="Yield_destroy")
+ use,intrinsic :: ISO_C_BINDING
+ implicit none
+ type(C_PTR), value :: this;
+end subroutine
+end interface
+end module Yield_I
diff --git a/sammy/src/salmon/interface/fortran/Yield_M.f90 b/sammy/src/salmon/interface/fortran/Yield_M.f90
new file mode 100644
index 0000000000000000000000000000000000000000..30eded3bbfee18acbe0a0afcafb62a2bf3732fa2
--- /dev/null
+++ b/sammy/src/salmon/interface/fortran/Yield_M.f90
@@ -0,0 +1,226 @@
+!>
+!! This file has been dynamically generated by Class Interface Xml (CIX)
+!! DO NOT MODIFY THIS FILE -- CHANGES WILL BE OVERWRITTEN UPON REGENERATION
+!! If changes need to occur, modify the appropriate CIX xml file
+!! Date Generated: Wed Apr 08 13:29:36 EDT 2020
+!! If any issues are experiences with this generated file that cannot be fixed
+!! with adjustment of the CIX xml file, please contact Robert A. Lefebvre, raq@ornl.gov
+!!/
+module Yield_M
+use, intrinsic :: ISO_C_BINDING
+use Yield_I
+type Yield
+ type(C_PTR) :: instance_ptr=C_NULL_PTR
+ contains
+ procedure, pass(this) :: setNormTypeSelfShielded => Yield_setNormTypeSelfShielded
+ procedure, pass(this) :: setNormTypeDivideByN => Yield_setNormTypeDivideByN
+ procedure, pass(this) :: setNormTypeTimesSigTot => Yield_setNormTypeTimesSigTot
+ procedure, pass(this) :: setCalcDerivs => Yield_setCalcDerivs
+ procedure, pass(this) :: setNumTotalParam => Yield_setNumTotalParam
+ procedure, pass(this) :: setY0 => Yield_setY0
+ procedure, pass(this) :: setY1 => Yield_setY1
+ procedure, pass(this) :: setY2 => Yield_setY2
+ procedure, pass(this) :: setTotY => Yield_setTotY
+ procedure, pass(this) :: setDerivY0 => Yield_setDerivY0
+ procedure, pass(this) :: setDerivY1 => Yield_setDerivY1
+ procedure, pass(this) :: setDerivY2 => Yield_setDerivY2
+ procedure, pass(this) :: setDerivTotY => Yield_setDerivTotY
+ procedure, pass(this) :: getNormString => Yield_getNormString
+ procedure, pass(this) :: normTypeIsSelfShielded => Yield_normTypeIsSelfShielded
+ procedure, pass(this) :: normTypeIsDivideByN => Yield_normTypeIsDivideByN
+ procedure, pass(this) :: normTypeIsTimesSigTot => Yield_normTypeIsTimesSigTot
+ procedure, pass(this) :: getCalcDerivs => Yield_getCalcDerivs
+ procedure, pass(this) :: getNumTotalParam => Yield_getNumTotalParam
+ procedure, pass(this) :: getY0 => Yield_getY0
+ procedure, pass(this) :: getY1 => Yield_getY1
+ procedure, pass(this) :: getY2 => Yield_getY2
+ procedure, pass(this) :: getTotY => Yield_getTotY
+ procedure, pass(this) :: getDerivY0 => Yield_getDerivY0
+ procedure, pass(this) :: getDerivY1 => Yield_getDerivY1
+ procedure, pass(this) :: getDerivY2 => Yield_getDerivY2
+ procedure, pass(this) :: getDerivTotY => Yield_getDerivTotY
+ procedure, pass(this) :: initialize => Yield_initialize
+ procedure, pass(this) :: destroy => Yield_destroy
+end type Yield
+contains
+subroutine Yield_setNormTypeSelfShielded(this)
+ implicit none
+ class(Yield)::this
+ call f_Yield_setNormTypeSelfShielded(this%instance_ptr)
+end subroutine
+subroutine Yield_setNormTypeDivideByN(this)
+ implicit none
+ class(Yield)::this
+ call f_Yield_setNormTypeDivideByN(this%instance_ptr)
+end subroutine
+subroutine Yield_setNormTypeTimesSigTot(this)
+ implicit none
+ class(Yield)::this
+ call f_Yield_setNormTypeTimesSigTot(this%instance_ptr)
+end subroutine
+subroutine Yield_setCalcDerivs(this, sb)
+ implicit none
+ class(Yield)::this
+ logical(C_BOOL)::sb
+ call f_Yield_setCalcDerivs(this%instance_ptr, sb)
+end subroutine
+subroutine Yield_setNumTotalParam(this, np)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::np
+ call f_Yield_setNumTotalParam(this%instance_ptr, np)
+end subroutine
+subroutine Yield_setY0(this, y0)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE)::y0
+ call f_Yield_setY0(this%instance_ptr, y0)
+end subroutine
+subroutine Yield_setY1(this, y1)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE)::y1
+ call f_Yield_setY1(this%instance_ptr, y1)
+end subroutine
+subroutine Yield_setY2(this, y2)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE)::y2
+ call f_Yield_setY2(this%instance_ptr, y2)
+end subroutine
+subroutine Yield_setTotY(this, toty)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE)::toty
+ call f_Yield_setTotY(this%instance_ptr, toty)
+end subroutine
+subroutine Yield_setDerivY0(this, iPar, dy0)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE)::dy0
+ call f_Yield_setDerivY0(this%instance_ptr, iPar-1,dy0)
+end subroutine
+subroutine Yield_setDerivY1(this, iPar, dy1)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE)::dy1
+ call f_Yield_setDerivY1(this%instance_ptr, iPar-1,dy1)
+end subroutine
+subroutine Yield_setDerivY2(this, iPar, dy2)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE)::dy2
+ call f_Yield_setDerivY2(this%instance_ptr, iPar-1,dy2)
+end subroutine
+subroutine Yield_setDerivTotY(this, iPar, dy2)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE)::dy2
+ call f_Yield_setDerivTotY(this%instance_ptr, iPar-1,dy2)
+end subroutine
+subroutine Yield_getNormString(this, string2Populate)
+ implicit none
+ class(Yield)::this
+ character(len=*) :: string2Populate
+ integer(C_INT) :: lengthOfString2Populate
+ lengthOfString2Populate = len(string2Populate)
+ call f_Yield_getNormString(this%instance_ptr, string2Populate, lengthOfString2Populate)
+end subroutine
+function Yield_normTypeIsSelfShielded(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ logical(C_BOOL):: result2Return
+ result2Return=f_Yield_normTypeIsSelfShielded(this%instance_ptr)
+end function
+function Yield_normTypeIsDivideByN(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ logical(C_BOOL):: result2Return
+ result2Return=f_Yield_normTypeIsDivideByN(this%instance_ptr)
+end function
+function Yield_normTypeIsTimesSigTot(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ logical(C_BOOL):: result2Return
+ result2Return=f_Yield_normTypeIsTimesSigTot(this%instance_ptr)
+end function
+function Yield_getCalcDerivs(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ logical(C_BOOL):: result2Return
+ result2Return=f_Yield_getCalcDerivs(this%instance_ptr)
+end function
+function Yield_getNumTotalParam(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ integer(C_INT):: result2Return
+ result2Return=f_Yield_getNumTotalParam(this%instance_ptr)
+end function
+function Yield_getY0(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getY0(this%instance_ptr)
+end function
+function Yield_getY1(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getY1(this%instance_ptr)
+end function
+function Yield_getY2(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getY2(this%instance_ptr)
+end function
+function Yield_getTotY(this) result(result2Return)
+ implicit none
+ class(Yield)::this
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getTotY(this%instance_ptr)
+end function
+function Yield_getDerivY0(this, iPar) result(result2Return)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getDerivY0(this%instance_ptr, iPar-1)
+end function
+function Yield_getDerivY1(this, iPar) result(result2Return)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getDerivY1(this%instance_ptr, iPar-1)
+end function
+function Yield_getDerivY2(this, iPar) result(result2Return)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getDerivY2(this%instance_ptr, iPar-1)
+end function
+function Yield_getDerivTotY(this, iPar) result(result2Return)
+ implicit none
+ class(Yield)::this
+ integer(C_INT)::iPar
+ real(C_DOUBLE):: result2Return
+ result2Return=f_Yield_getDerivTotY(this%instance_ptr, iPar-1)
+end function
+subroutine Yield_initialize(this)
+ implicit none
+ class(Yield) :: this
+ this%instance_ptr = f_Yield_initialize()
+end subroutine
+subroutine Yield_destroy(this)
+ implicit none
+ class(Yield) :: this
+ call f_Yield_destroy(this%instance_ptr)
+ this%instance_ptr = C_NULL_PTR
+end subroutine
+end module Yield_M
diff --git a/sammy/src/salmon/tests/CMakeLists.txt b/sammy/src/salmon/tests/CMakeLists.txt
index 7f3bf5f575aa5a84ae422559833f4daba26b95e1..9bbab0a950d21215703c8635e40f553096218207 100644
--- a/sammy/src/salmon/tests/CMakeLists.txt
+++ b/sammy/src/salmon/tests/CMakeLists.txt
@@ -2,3 +2,4 @@ include(NemesisTest)
ADD_NEMESIS_TEST(ExperimentalParametersTest.cpp TIMEOUT 60000 NP 1 )
ADD_NEMESIS_TEST(MultipleScatteringTest.cpp TIMEOUT 60000 NP 1 )
+ADD_NEMESIS_TEST(YieldTest.cpp TIMEOUT 60000 NP 1 )
diff --git a/sammy/src/salmon/tests/ExperimentalParametersTest.cpp b/sammy/src/salmon/tests/ExperimentalParametersTest.cpp
index 2d6d9c19232b4c575eb32a6fef82f7757428b875..96f48c0b2359ad81840dd7acb95e3def85326f64 100644
--- a/sammy/src/salmon/tests/ExperimentalParametersTest.cpp
+++ b/sammy/src/salmon/tests/ExperimentalParametersTest.cpp
@@ -11,20 +11,23 @@ TEST(SampleDimensions,setterGetter){
sammy::SampleDimensions result;
- result.setRadius( 1.0);
- result.setHeight( 2.0);
- result.setThickness( 3.0);
- result.setBeamRadius(4.0);
-
- double goldRadius = 1.0;
- double goldHeight = 2.0;
- double goldThickness = 3.0;
- double goldBeamRadius = 4.0;
-
- ASSERT_EQ(goldRadius,result.getRadius());
- ASSERT_EQ(goldHeight,result.getHeight());
- ASSERT_EQ(goldThickness,result.getThickness());
- ASSERT_EQ(goldBeamRadius,result.getBeamRadius());
+ result.setRadius( 1.0);
+ result.setHeight( 2.0);
+ result.setThickness( 3.0);
+ result.setArealDensity(4.0);
+ result.setBeamRadius( 5.0);
+
+ double goldRadius = 1.0;
+ double goldHeight = 2.0;
+ double goldThickness = 3.0;
+ double goldArealDensity = 4.0;
+ double goldBeamRadius = 5.0;
+
+ ASSERT_EQ(goldRadius, result.getRadius());
+ ASSERT_EQ(goldHeight, result.getHeight());
+ ASSERT_EQ(goldThickness, result.getThickness());
+ ASSERT_EQ(goldArealDensity,result.getArealDensity());
+ ASSERT_EQ(goldBeamRadius, result.getBeamRadius());
}
/**********************************************************
diff --git a/sammy/src/salmon/tests/YieldTest.cpp b/sammy/src/salmon/tests/YieldTest.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..3ce0bb696b0cb42f192899f21150fdc9f3662354
--- /dev/null
+++ b/sammy/src/salmon/tests/YieldTest.cpp
@@ -0,0 +1,58 @@
+#include "Nemesis/gtest/Gtest_Functions.hh"
+#include "Nemesis/gtest/nemesis_gtest.hh"
+#include "../Yield.h"
+
+
+
+/**************************************************************************************************
+* Test the getting and setting functions
+**************************************************************************************************/
+TEST(Yield,setterGetter){
+
+ sammy::Yield yield;
+
+ int gold_def_normType_NOTSET = true;
+ int gold_normType_tst = true; // <-- Normtype is now set to normalization::TIMESSIGTOT
+ double gold_Y0 = 0.1;
+ double gold_Y1 = 0.2;
+ double gold_Y2 = 0.3;
+ double gold_derivY0 = 0.1;
+ double gold_derivY1 = 0.2;
+ double gold_derivY2 = 0.3;
+ double gold_totY = gold_Y0 + gold_Y1 + gold_Y2;
+ double gold_derivTotY = -0.1;
+ bool gold_def_calcDerivs = false;
+ bool gold_calcDerivs = true;
+ int gold_numParam = 20;
+
+
+ yield.setY0( gold_Y0);
+ yield.setY1( gold_Y1);
+ yield.setY2( gold_Y2);
+ yield.setTotY( gold_totY);
+ yield.setDerivY0(0, gold_derivY0);
+ yield.setDerivY1(0, gold_derivY1);
+ yield.setDerivY2(0, gold_derivY2);
+ yield.setDerivTotY(0, gold_derivTotY);
+ yield.setNumTotalParam( gold_numParam);
+
+
+ ASSERT_EQ(gold_def_normType_NOTSET, yield.normTypeIsNotSet()); // <-- Default normalization type
+ ASSERT_EQ(gold_def_calcDerivs, yield.getCalcDerivs());
+
+ yield.setNormTypeTimesSigTot();
+ yield.setCalcDerivs(gold_calcDerivs);
+
+ ASSERT_EQ(gold_normType_tst, yield.normTypeIsTimesSigTot());
+ ASSERT_EQ(gold_Y0, yield.getY0());
+ ASSERT_EQ(gold_Y1, yield.getY1());
+ ASSERT_EQ(gold_Y2, yield.getY2());
+ ASSERT_EQ(gold_derivY0, yield.getDerivY0(0));
+ ASSERT_EQ(gold_derivY1, yield.getDerivY1(0));
+ ASSERT_EQ(gold_derivY2, yield.getDerivY2(0));
+ ASSERT_EQ(gold_totY, yield.getTotY());
+ ASSERT_EQ(gold_derivTotY, yield.getDerivTotY(0));
+ ASSERT_EQ(gold_numParam, yield.getNumTotalParam());
+ ASSERT_EQ(gold_calcDerivs, yield.getCalcDerivs());
+
+}
diff --git a/sammy/src/sam/msam.F b/sammy/src/sam/msam.F
index 602134bd4721930bc2018ab44aeeb35dda212ecf..cce141037281bc60565fcd220c8d4bb0e27f9b01 100755
--- a/sammy/src/sam/msam.F
+++ b/sammy/src/sam/msam.F
@@ -16,6 +16,7 @@ C
use EndfData_common_m
use MultScatPars_common_m
use ExpPars_common_m
+ use Observable_common_m
use ssm_m
use Samdat_0_M
use Samodf_0_m
@@ -43,6 +44,7 @@ C
call multScat%initialize()
call sampleDim%initialize()
call covData%initialize()
+ call yld%initialize()
call covData%makeNewUCovariance()
call covData%makeNewCovariance()
call covData%makeNewInverseUCovariance()
diff --git a/sammy/src/sammy/CMakeLists.txt b/sammy/src/sammy/CMakeLists.txt
index 8eae11683ba0c47d0e0dbdccbeddd82a79f23f83..bc1bfa2dd475703475dc1e85d6e5a1887ec259a3 100644
--- a/sammy/src/sammy/CMakeLists.txt
+++ b/sammy/src/sammy/CMakeLists.txt
@@ -572,7 +572,8 @@ APPEND_SET(SAMMY_SOURCES
../blk/AllocateFunctions.f90
../blk/ExpPars_common.f90
../blk/MultScatPars_common.f90
- ../blk/CapYCorrections_common.f90
+ ../blk/Convolution_common.f90
+ ../blk/Observable_common.f90
)
TRIBITS_ADD_LIBRARY(
diff --git a/sammy/src/ssm/mssm04.f90 b/sammy/src/ssm/mssm04.f90
index 12a4ea9c1839210d790895cf1f15c0dfcf21e615..df05f554aa532d8b8d2c47a86f7e0b556f5062ea 100644
--- a/sammy/src/ssm/mssm04.f90
+++ b/sammy/src/ssm/mssm04.f90
@@ -40,6 +40,7 @@ module ssm_4_m
Itimes = 1
Kountr = 0
!
+! JMB: The a and b here are Y's at energy a and b read from a file
IF (Mcy2.EQ.1) CALL Gety2zzz (Emy2_a, Emy2_b, Y2tab_a, Y2tab_b)
!
CALL Findpr (multScat%getNumTheta(), Knthet)
diff --git a/sammy/src/ssm/mssm19.f90 b/sammy/src/ssm/mssm19.f90
index 7449e786a619c416cce49db50e05dde70d2ac8a9..f885014f8df3e9d1e66e03a338f00e992eb0c54f 100644
--- a/sammy/src/ssm/mssm19.f90
+++ b/sammy/src/ssm/mssm19.f90
@@ -11,32 +11,47 @@ module ssm_19_m
use ifwrit_m
use fixedr_m
use logic_ssm_common_m
- use CapYCorrections_common_m
+ use Convolution_common_m
+ use Observable_common_m
+ use MultScatPars_common_m
+ use ExpPars_common_m
+ use Yield_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
+ logical(C_BOOL) calc_derivs
DIMENSION Capsig(*), Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Sigxxx(*), Dasigx(*), Dbsigx(*)
DATA Zero /0.0d0/, One /1.0d0/
call capYieldCor%initialize()
+ call yNorm%initialize()
!
! *** Default normalization is such that Y0 = (1-e)Capture /(n Total)
-! *** which corresponds to Kyield=1 ("cross section", not "yield")
+! *** which corresponds to yld%normTypeIsDivideByN()==true ("cross section", not "yield")
!
IF (Noyzer.NE.0) Y0 = Zero
+
+ ! -----------------------------------
+ ! Put things into C++ memory here
+ ! -----------------------------------
+ call yld%setY0(Y0)
+ calc_derivs = Ksolve.NE.2
+ call yld%setCalcDerivs(calc_derivs)
+ call multScat%setNumScatters(0)
+ call yld%setNumTotalParam(Nnpar)
+ call sampleDim%setArealDensity(Dthick)
!
Cap = Capsig(Nn)
- IF (Kyield.EQ.0) THEN
- Y0 = Y0*Dthick
- ELSE IF (Kyield.EQ.1) THEN
- ELSE IF (Kyield.EQ.2) THEN
- Y0 = Y0*Total*Dthick
- END IF
+
+ ! -----------------------------------
+ ! Now normalize the yield
+ ! -----------------------------------
+ call yNorm%norm(yld,sampleDim,multScat,Total)
!
- Sigxxx(1) = Y0
+ call yld%setTotY(yld%getY0())
+ Sigxxx(1) = yld%getTotY()
+ Y0 = yld%getY0()
!
! If the gamma attenuation factor C1 is not zero
IF (Gam_att_C1.NE.0.0) THEN
- dC1 = dGam_att_C1
- Td = Total*Dthick
! Sigxxx(1) is the place we store the output to ODF/LST files
! -- Here we are modifying Sigxxx(1) with the gamma attenuation
! -- correction
@@ -61,13 +76,13 @@ module ssm_19_m
IF (Kvthck.GT.0) THEN
! *** Add deriv of Y0 wrt Thickness
Kv = Kvthck - Nvadif - Ndasig
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
IF (Noyzer.EQ.0) THEN
Dbsigx(Kv) = Exp1*Cap + Dbsigx(Kv)
! ELSE
! disabling Y0 portion of calculation
END IF
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Sn = Dthick*Total
IF (Noyzer.EQ.0) THEN
IF (Sn.LT.0.3d0) THEN
@@ -81,18 +96,21 @@ module ssm_19_m
F = Zero
END IF
Dbsigx(Kv) = F + Dbsigx(Kv)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
IF (Noyzer.EQ.0) Dbsigx(Kv)= Exp1*Cap*Total + Dbsigx(Kv)
END IF
END IF
!
IF (Sensin.NE.Zero .AND. Ksensc.GT.0) THEN
- IF (Kyield.EQ.0) Asensn = Asensn * Dthick
- IF (Kyield.EQ.2) Asensn = Asensn * Dthick * Total
+ IF (yld%normTypeIsSelfShielded()) Asensn = Asensn * Dthick
+ IF (yld%normTypeIsTimesSigTot()) Asensn = Asensn * Dthick * Total
! *** Add deriv of Y0 wrt neutron sensitivity multiplier
Kv = Ksensc - Nvadif - Ndasig
Dbsigx(Kv) = Asensn
END IF
+
+ call capYieldCor%destroy()
+
RETURN
END
!
@@ -104,7 +122,8 @@ module ssm_19_m
use fixedi_m
use ifwrit_m
use logic_ssm_common_m
- use CapYCorrections_common_m
+ use Convolution_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
DIMENSION Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Ds(*)
DATA One /1.0d0/
@@ -120,18 +139,18 @@ module ssm_19_m
IF (Noyzer.EQ.0) THEN
! *** First, the derivative of Y0 wrt u(Iipar)
! *** (Nn=Energy) JMB
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total * Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
IF (Iipar.EQ.1) Dddyy0 = Dddyy0 + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + (One-Exp1)*Dcapsi(Iipar,Nn) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap*Dthick
END IF
diff --git a/sammy/src/ssm/mssm20.f90 b/sammy/src/ssm/mssm20.f90
index 8ce787e297afb0aabad5e1ca0c9610bdeffbae10..9c3c3ce709cba51b2e659d7cdd208d7f2c60f783 100644
--- a/sammy/src/ssm/mssm20.f90
+++ b/sammy/src/ssm/mssm20.f90
@@ -149,54 +149,72 @@ module ssm_20_m
use ifwrit_m
use fixedr_m
use logic_ssm_common_m
- use CapYCorrections_common_m
+ use Convolution_common_m
+ use Observable_common_m
+ use MultScatPars_common_m
+ use ExpPars_common_m
+ use Yield_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
+ logical(C_BOOL) calc_derivs
DIMENSION Capsig(*), Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*), &
Sigxxx(*), Dasigx(*), Dbsigx(*)
DATA Zero /0.0d0/, One /1.0d0/
call capYieldCor%initialize()
+ call yNorm%initialize()
!
! *** Default normalization is such that Y0 = (1-e)Capture /(n Total)
-! *** which corresponds to Kyield=1 ("cross section", not "yield")
+! *** which corresponds to yld%normTypeIsDivideByN()==true ("cross section", not "yield")
!
IF (Noyzer.NE.0) Y0 = Zero
+
+ ! -----------------------------------
+ ! Put things into C++ memory here
+ ! -----------------------------------
+ call yld%setY0(Y0)
+ call yld%setY1(Y1)
+ if( Ksolve.ne.2 ) then
+ do j=1,Nnpar
+ call yld%setDerivY1(J,Dy1(J))
+ end do
+ end if
+ calc_derivs = Ksolve.NE.2
+ call yld%setCalcDerivs(calc_derivs)
+ call multScat%setNumScatters(1)
+ call yld%setNumTotalParam(Nnpar)
+ call sampleDim%setArealDensity(Dthick)
+
+ ! -----------------------------------
+ ! Now normalize the yield
+ ! -----------------------------------
+ call yNorm%norm(yld,sampleDim,multScat,Total)
!
Cap = Capsig(Nn)
- IF (Kyield.EQ.0) THEN
- Y0 = Y0*Dthick
- ELSE IF (Kyield.EQ.1) THEN
- Y1 = Y1/Dthick
- IF (Ksolve.NE.2) THEN
- DO J=1,Nnpar
- Dy1(J) = Dy1(J)/Dthick
- END DO
- END IF
- ELSE IF (Kyield.EQ.2) THEN
- Y0 = Y0*Total*Dthick
- Y1 = Y1*Total
- IF (Ksolve.NE.2) THEN
- DO J=1,Nnpar
- Dy1(J) = Dy1(J)*Total
- END DO
- END IF
- END IF
-!
-!
- Sigxxx(1) = Y0
- Sigxxx(1) = Sigxxx(1) + Y1
+
+ ! -----------------------------------
+ ! Return the values to the old places
+ ! -----------------------------------
+ call yld%setTotY(yld%getY0() + yld%getY1())
+ Sigxxx(1) = yld%getTotY()
+ Y0 = yld%getY0()
+ Y1 = yld%getY1()
+ if( yld%getCalcDerivs() ) then
+ do j=1,yld%getNumTotalParam()
+ Dy1(j) = yld%getDerivY1(j)
+ end do
+ end if
!
! If the gamma attenuation factor C1 is not zero
IF (Gam_att_C1.NE.0.0) THEN
- dC1 = dGam_att_C1
- Td = Total*Dthick
! Sigxxx(1) is the place we store the output to ODF/LST files
! -- Here we are modifying Sigxxx(1) with the gamma attenuation
! -- correction
call capYieldCor%gammaAttenCorr( Dthick,Total,Gam_att_C1,Sigxxx(1) )
END IF
!
-!
+ ! -----------------------------------
+ ! Write yields to file
+ ! -----------------------------------
IF (Kssmpr.EQ.1) THEN
WRITE (60,12000) Em, Y0, Y1, Sigxxx(1)
12000 FORMAT (1X, f15.8, 1P5E15.6)
@@ -213,13 +231,13 @@ module ssm_20_m
IF (Kvthck.GT.0) THEN
! *** Add deriv of Y0 wrt Thickness; modify deriv of Y1 & Y2 wrt Thick
Kv = Kvthck - Nvadif - Ndasig
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
IF (Noyzer.EQ.0) THEN
Dbsigx(Kv) = Exp1*Cap + Dbsigx(Kv)
! ELSE
! disabling Y0 portion of calculation
END IF
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Sn = Dthick*Total
IF (Noyzer.EQ.0) THEN
IF (Sn.LT.0.3d0) THEN
@@ -233,18 +251,21 @@ module ssm_20_m
F = Zero
END IF
Dbsigx(Kv) = F + Dbsigx(Kv) - Y1/Dthick
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
IF (Noyzer.EQ.0) Dbsigx(Kv)= Exp1*Cap*Total + Dbsigx(Kv)
END IF
END IF
!
IF (Sensin.NE.Zero .AND. Ksensc.GT.0) THEN
- IF (Kyield.EQ.0) Asensn = Asensn * Dthick
- IF (Kyield.EQ.2) Asensn = Asensn * Dthick * Total
+ IF (yld%normTypeIsSelfShielded()) Asensn = Asensn * Dthick
+ IF (yld%normTypeIsTimesSigTot()) Asensn = Asensn * Dthick * Total
! *** Add deriv of Y0 wrt neutron sensitivity multiplier
Kv = Ksensc - Nvadif - Ndasig
Dbsigx(Kv) = Asensn
END IF
+
+ call capYieldCor%destroy()
+
RETURN
END
!
@@ -256,7 +277,8 @@ module ssm_20_m
use fixedi_m
use ifwrit_m
use logic_ssm_common_m
- use CapYCorrections_common_m
+ use Convolution_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
DIMENSION Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*), Ds(*)
DATA One /1.0d0/
@@ -272,29 +294,29 @@ module ssm_20_m
IF (Noyzer.EQ.0) THEN
! *** First, the derivative of Y0 wrt u(Iipar)
! *** (Nn=Energy) JMB
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total * Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
IF (Iipar.EQ.1) Dddyy0 = Dddyy0 + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + (One-Exp1)*Dcapsi(Iipar,Nn) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap*Dthick
END IF
END IF
!
! *** Next, the derivative of Y1 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar) + &
Y1*Dtotsi(Iipar,Nn)/Total
END IF
diff --git a/sammy/src/ssm/mssm21.f90 b/sammy/src/ssm/mssm21.f90
index 49047ddebc76f8384bd54c303d5ab3eeabacda7d..0bcd285b693330649a16d1b42a3f57700a099279 100644
--- a/sammy/src/ssm/mssm21.f90
+++ b/sammy/src/ssm/mssm21.f90
@@ -49,9 +49,16 @@ module ssm_21_m
use fixedi_m
use ifwrit_m
use fixedr_m
- use CapYCorrections_common_m
+ use Convolution_common_m
use logic_ssm_common_m
+ use Observable_common_m
+ use MultScatPars_common_m
+ use ExpPars_common_m
+ use EndfData_common_m
+ use Yield_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
+ logical(C_BOOL) calc_derivs
+ double precision, dimension(:), allocatable :: Dabsigx
DIMENSION Capsig(*), Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*),&
Dy2(*), Sigxxx(*), Dasigx(*), Dbsigx(*)
DATA Zero /0.0d0/, One /1.0d0/
@@ -59,54 +66,72 @@ module ssm_21_m
!type(CapYieldCorrections)::capYieldCor
call capYieldCor%initialize()
+ call yNorm%initialize()
+
+ ! ------- Nnpar, Ndasig, and Ndbsig are set to zero when Ksolve = 2 -----------
+
!
! *** Default normalization is such that Y0 = (1-e)Capture /(n Total)
-! *** which corresponds to Kyield=1 ("cross section", not "yield")
+! *** which corresponds to yld%normTypeIsDivideByN()==true ("cross section", not "yield")
!
IF (Noyzer.NE.0) Y0 = Zero
+
+ ! -----------------------------------
+ ! Put things into C++ memory here
+ ! -----------------------------------
+ call yld%setY0(Y0)
+ call yld%setY1(Y1)
+ call yld%setY2(Y2)
+ if( Ksolve.ne.2 ) then
+ do j=1,Nnpar
+ call yld%setDerivY1(J,Dy1(J))
+ call yld%setDerivY2(J,Dy2(J))
+ end do
+ end if
+ calc_derivs = Ksolve.NE.2
+ call yld%setCalcDerivs(calc_derivs)
+ call multScat%setNumScatters(2)
+ call yld%setNumTotalParam(Nnpar)
+ call sampleDim%setArealDensity(Dthick)
+
+ ! -----------------------------------
+ ! Now normalize the yield
+ ! -----------------------------------
+ call yNorm%norm(yld,sampleDim,multScat,Total)
!
Cap = Capsig(Nn)
- IF (Kyield.EQ.0) THEN
- Y0 = Y0*Dthick
- ELSE IF (Kyield.EQ.1) THEN
- Y1 = Y1/Dthick
- Y2 = Y2/Dthick
- IF (Ksolve.NE.2) THEN
- DO J=1,Nnpar
- Dy1(J) = Dy1(J)/Dthick
- END DO
- DO J=1,Nnpar
- Dy2(J) = Dy2(J)/Dthick
- END DO
- END IF
- ELSE IF (Kyield.EQ.2) THEN
- Y0 = Y0*Total*Dthick
- Y1 = Y1*Total
- Y2 = Y2*Total
- IF (Ksolve.NE.2) THEN
- DO J=1,Nnpar
- Dy1(J) = Dy1(J)*Total
- END DO
- DO J=1,Nnpar
- Dy2(J) = Dy2(J)*Total
- END DO
- END IF
- END IF
- IF (Mcy2.NE.0) Y2 = Y2tab
-!
- Sigxxx(1) = Y0 + Y1 + Y2
+
+ IF (Mcy2.NE.0) then
+ Y2 = Y2tab
+ call yld%setY2(Y2tab)
+ end if
+!
+ ! -----------------------------------
+ ! Return the values to the old places
+ ! -----------------------------------
+ call yld%setTotY(yld%getY0() + yld%getY1() + yld%getY2())
+ Sigxxx(1) = yld%getTotY()
+ Y0 = yld%getY0()
+ Y1 = yld%getY1()
+ Y2 = yld%getY2()
+ if( yld%getCalcDerivs() ) then
+ do j=1,yld%getNumTotalParam()
+ Dy1(j) = yld%getDerivY1(j)
+ Dy2(j) = yld%getDerivY2(j)
+ end do
+ end if
!
! If the gamma attenuation factor C1 is not zero
IF (Gam_att_C1.NE.0.0) THEN
- dC1 = dGam_att_C1
- Td = Total*Dthick
! Sigxxx(1) is the place we store the output to ODF/LST files
! -- Here we are modifying Sigxxx(1) with the gamma attenuation
! -- correction
call capYieldCor%gammaAttenCorr( Dthick,Total,Gam_att_C1,Sigxxx(1) )
END IF
!
-!
+ ! -----------------------------------
+ ! Write yields to file
+ ! -----------------------------------
IF (Kssmpr.EQ.1) THEN
WRITE (60,12000) Em, Y0, Y1, Y2, Sigxxx(1)
12000 FORMAT (1X, f15.8, 1P5E15.6)
@@ -115,7 +140,7 @@ module ssm_21_m
IF (Ksolve.EQ.2) RETURN
! ######################################
! It seems that Ndasig is a subset of varied parameters and that
-! Ndbsig is the total number of varied parameters minus Ndasig
+! Ndbsig is another subset of varied parameters minus Ndasig
! ######################################
IF (Ndasig.GT.0) CALL Ynrm_2_Derivative (Cap, Dtotsi, Dcapsi, &
Dasigx, Dy1, Dy2, Y0, Y1, Y2, Dthick, Total, Exp1, Dddyy0, &
@@ -127,13 +152,13 @@ module ssm_21_m
IF (Kvthck.GT.0) THEN
! *** Add deriv of Y0 wrt Thickness; modify deriv of Y1 & Y2 wrt Thick
Kv = Kvthck - Nvadif - Ndasig
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
IF (Noyzer.EQ.0) THEN
Dbsigx(Kv) = Exp1*Cap + Dbsigx(Kv)
! ELSE
! disabling Y0 portion of calculation
END IF
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Sn = Dthick*Total
IF (Noyzer.EQ.0) THEN
IF (Sn.LT.0.3d0) THEN
@@ -148,18 +173,21 @@ module ssm_21_m
END IF
Dbsigx(Kv) = F + Dbsigx(Kv) - Y1/Dthick
IF (Double_Plus) Dbsigx(Kv) = Dbsigx(Kv)-Y2/Dthick
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
IF (Noyzer.EQ.0) Dbsigx(Kv)= Exp1*Cap*Total + Dbsigx(Kv)
END IF
END IF
!
IF (Sensin.NE.Zero .AND. Ksensc.GT.0) THEN
- IF (Kyield.EQ.0) Asensn = Asensn * Dthick
- IF (Kyield.EQ.2) Asensn = Asensn * Dthick * Total
+ IF (yld%normTypeIsSelfShielded()) Asensn = Asensn * Dthick
+ IF (yld%normTypeIsTimesSigTot()) Asensn = Asensn * Dthick * Total
! *** Add deriv of Y0 wrt neutron sensitivity multiplier
Kv = Ksensc - Nvadif - Ndasig
Dbsigx(Kv) = Asensn
END IF
+
+ call capYieldCor%destroy()
+
RETURN
END
!
@@ -171,10 +199,12 @@ module ssm_21_m
use fixedi_m
use ifwrit_m
use logic_ssm_common_m
- use CapYCorrections_common_m
+ use Convolution_common_m
+ use Observable_common_m
IMPLICIT DOUBLE PRECISION (a-h,o-z)
DIMENSION Dtotsi(Nnparx,*), Dcapsi(Nnparx,*), Dy1(*), Dy2(*),Ds(*)
DATA One /1.0d0/
+
!
! Dthick = N [at/b]
! Cap = sigma_gamma
@@ -187,39 +217,39 @@ module ssm_21_m
IF (Noyzer.EQ.0) THEN
! *** First, the derivative of Y0 wrt u(Iipar)
! *** (Nn=Energy) JMB
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total * Dthick
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
IF (Iipar.EQ.1) Dddyy0 = Dddyy0 + &
Y0* (Dcapsi(Iipar,Nn)/Cap-Dtotsi(Iipar,Nn)/Total) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap/Total
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + (One-Exp1)*Dcapsi(Iipar,Nn) &
+ Dtotsi(Iipar,Nn)*Exp1*Cap*Dthick
END IF
END IF
!
! *** Next, the derivative of Y1 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + Dy1(Iipar) + &
Y1*Dtotsi(Iipar,Nn)/Total
END IF
!
! *** Finally, the derivative of Y2 wrt u(Iipar)
- IF (Kyield.EQ.0) THEN
+ IF (yld%normTypeIsSelfShielded()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar)
- ELSE IF (Kyield.EQ.1) THEN
+ ELSE IF (yld%normTypeIsDivideByN()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar)
- ELSE IF (Kyield.EQ.2) THEN
+ ELSE IF (yld%normTypeIsTimesSigTot()) THEN
Ds(Ipar) = Ds(Ipar) + Dy2(Iipar) &
+ Y2*Dtotsi(Iipar,Nn)/Total
END IF