diff --git a/src/endf.cpp b/src/endf.cpp index 859183d12..8a63893a7 100644 --- a/src/endf.cpp +++ b/src/endf.cpp @@ -144,4 +144,35 @@ double Tabulated1D::operator()(double x) const } } +//============================================================================== +// CoherentElasticXS implementation +//============================================================================== + +CoherentElasticXS::CoherentElasticXS(hid_t dset) +{ + // Read 2D array from dataset + xt::xarray arr; + read_dataset(dset, arr); + + // Get views for Bragg edges and structure factors + auto E = xt::view(arr, 0); + auto s = xt::view(arr, 1); + + // Copy Bragg edges and partial sums of structure factors + std::copy(E.begin(), E.end(), std::back_inserter(bragg_edges_)); + std::copy(s.begin(), s.end(), std::back_inserter(factors_)); +} + +double CoherentElasticXS::operator()(double E) const +{ + if (E < bragg_edges_[0]) { + // If energy is below that of the lowest Bragg peak, the elastic cross + // section will be zero + return 0.0; + } else { + auto i_grid = lower_bound_index(bragg_edges_.begin(), bragg_edges_.end(), E); + return factors_[i_grid] / E; + } +} + } // namespace openmc diff --git a/src/endf.h b/src/endf.h index d2d67685d..369f5d5a6 100644 --- a/src/endf.h +++ b/src/endf.h @@ -73,6 +73,19 @@ private: std::vector y_; //!< values of ordinate }; +//============================================================================== +//! Coherent elastic scattering data from a crystalline material +//============================================================================== + +class CoherentElasticXS : public Function1D { + explicit CoherentElasticXS(hid_t dset); + double operator()(double E) const; +private: + std::vector bragg_edges_; //!< Bragg edges in [eV] + std::vector factors_; //!< Partial sums of structure factors [eV-b] +}; + + } // namespace openmc #endif // OPENMC_ENDF_H