Mantid
Loading...
Searching...
No Matches
ThermalNeutronBk2BkExpAlpha.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
10
11#include <cmath>
12#include <gsl/gsl_sf_erf.h>
13
14using namespace std;
15
16using namespace Mantid;
17
18using namespace Mantid::API;
19
20using namespace Mantid::CurveFitting;
21
23
24using namespace CurveFitting;
25
27
28//----------------------------------------------------------------------------------------------
32 // Geometry related
33 declareParameter("Width", 1.0);
34 declareParameter("Tcross", 1.0);
35
36 declareParameter("Alph0", 0.0);
37 declareParameter("Alph1", 0.0);
38 declareParameter("Alph0t", 0.0);
39 declareParameter("Alph1t", 0.0);
40}
41
42//----------------------------------------------------------------------------------------------
45void ThermalNeutronBk2BkExpAlpha::function1D(double *out, const double *xValues, const size_t nData) const {
46 double width = getParameter("Width");
47 double tcross = getParameter("Tcross");
48 double alph0 = getParameter("Alph0");
49 double alph1 = getParameter("Alph1");
50 double alph0t = getParameter("Alph0t");
51 double alph1t = getParameter("Alph1t");
52
53 for (size_t i = 0; i < nData; ++i) {
54 out[i] = corefunction(xValues[i], width, tcross, alph0, alph1, alph0t, alph1t);
55 }
56}
57
61 calNumericalDeriv(domain, jacobian);
62}
63
64//----------------------------------------------------------------------------------------------
67double ThermalNeutronBk2BkExpAlpha::corefunction(double dh, double width, double tcross, double alph0, double alph1,
68 double alph0t, double alph1t) const {
69 double n = 0.5 * gsl_sf_erfc(width * (tcross - 1.0 / dh));
70 double alpha = 1.0 / (n * (alph0 + alph1 * dh) + (1.0 - n) * (alph0t - alph1t / dh));
71
72 return alpha;
73}
74
75} // namespace Mantid::CurveFitting::Functions
#define DECLARE_FUNCTION(classname)
Macro for declaring a new type of function to be used with the FunctionFactory.
Base class that represents the domain of a function.
void calNumericalDeriv(const FunctionDomain &domain, Jacobian &jacobian)
Calculate numerical derivatives.
Definition: IFunction.cpp:1031
Represents the Jacobian in IFitFunction::functionDeriv.
Definition: Jacobian.h:22
double getParameter(size_t i) const override
Get i-th parameter.
ThermalNeutronBk2BkExpAlpha : Function to calculate Alpha of Bk2Bk Exponential function from Thermal ...
void function1D(double *out, const double *xValues, const size_t nData) const override
Override.
double corefunction(double dh, double width, double tcross, double alph0, double alph1, double alph0t, double alph1t) const
Core function (inline) to calcualte TOF_h from d-spacing.
void functionDeriv(const API::FunctionDomain &domain, API::Jacobian &jacobian) override
Derivative to overwrite.
Helper class which provides the Collimation Length for SANS instruments.
STL namespace.