Mantid
Loading...
Searching...
No Matches
ThermalNeutronBk2BkExpBeta.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::CurveFitting;
19
20using namespace Mantid::API;
21
23
24using namespace CurveFitting;
25
27
28//----------------------------------------------------------------------------------------------
32 declareParameter("Width", 1.0);
33 declareParameter("Tcross", 1.0);
34
35 declareParameter("Beta0", 0.0);
36 declareParameter("Beta1", 0.0);
37 declareParameter("Beta0t", 0.0);
38 declareParameter("Beta1t", 0.0);
39}
40
41//----------------------------------------------------------------------------------------------
44void ThermalNeutronBk2BkExpBeta::function1D(double *out, const double *xValues, const size_t nData) const {
45 double width = getParameter("Width");
46 double tcross = getParameter("Tcross");
47 double beta0 = getParameter("Beta0");
48 double beta1 = getParameter("Beta1");
49 double beta0t = getParameter("Beta0t");
50 double beta1t = getParameter("Beta1t");
51
52 for (size_t i = 0; i < nData; ++i) {
53 out[i] = corefunction(xValues[i], width, tcross, beta0, beta1, beta0t, beta1t);
54 }
55}
56
60 calNumericalDeriv(domain, jacobian);
61}
62
63//----------------------------------------------------------------------------------------------
66double ThermalNeutronBk2BkExpBeta::corefunction(double dh, double width, double tcross, double beta0, double beta1,
67 double beta0t, double beta1t) const {
68 double n = 0.5 * gsl_sf_erfc(width * (tcross - 1.0 / dh));
69 double beta = 1.0 / (n * (beta0 + beta1 * dh) + (1.0 - n) * (beta0t - beta1t / dh));
70
71 return beta;
72}
73
74} // 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.
ThermalNeutronBk2BkExpBETA : Function to calculate Beta of Bk2Bk Exponential function from Thermal Ne...
void function1D(double *out, const double *xValues, const size_t nData) const override
Override.
double corefunction(double dh, double width, double tcross, double beta0, double beta1, double beta0t, double beta1t) 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.