Skip to content

Commit

Permalink
addition of Buckingham potential
Browse files Browse the repository at this point in the history
  • Loading branch information
@Danbr4d
Danbr4d committed Apr 25, 2024
1 parent c8a1e46 commit 408996d
Show file tree
Hide file tree
Showing 4 changed files with 75 additions and 38 deletions.
6 changes: 4 additions & 2 deletions src/classes/shortRangeFunctions.cpp
View file
Expand Up @@ -9,7 +9,8 @@ EnumOptions<ShortRangeFunctions::Form> ShortRangeFunctions::forms()
return EnumOptions<ShortRangeFunctions::Form>("ShortRangeType",
{{ShortRangeFunctions::Form::None, "None"},
{ShortRangeFunctions::Form::LennardJones, "LJ", 2, 2},
{ShortRangeFunctions::Form::LennardJonesGeometric, "LJGeometric", 2, 2}});
{ShortRangeFunctions::Form::LennardJonesGeometric, "LJGeometric", 2, 2},
{ShortRangeFunctions::Form::Buckingham, "Buckingham", 3, 3}});
}

// Return parameters for specified form
Expand All @@ -18,7 +19,8 @@ const std::vector<std::string> &ShortRangeFunctions::parameters(Form form)
static std::map<ShortRangeFunctions::Form, std::vector<std::string>> params_ = {
{ShortRangeFunctions::Form::None, {}},
{ShortRangeFunctions::Form::LennardJones, {"epsilon", "sigma"}},
{ShortRangeFunctions::Form::LennardJonesGeometric, {"epsilon", "sigma"}}};
{ShortRangeFunctions::Form::LennardJonesGeometric, {"epsilon", "sigma"}},
{ShortRangeFunctions::Form::Buckingham, {"A", "B", "C"}}};
return params_[form];
}

Expand Down
7 changes: 4 additions & 3 deletions src/classes/shortRangeFunctions.h
View file
Expand Up @@ -14,9 +14,10 @@ class ShortRangeFunctions
public:
enum class Form
{
None, /* No short-range dispersive forces */
LennardJones, /* Lennard-Jones 12-6 form with Lorentz-Berthelot combination rules */
LennardJonesGeometric /* Lennard-Jones 12-6 form with Geometric combination rules */
None, /* No short-range dispersive forces */
LennardJones, /* Lennard-Jones 12-6 form with Lorentz-Berthelot combination rules */
LennardJonesGeometric, /* Lennard-Jones 12-6 form with Geometric combination rules */
Buckingham /* Buckingham form */
};
// Return enum options for form
static EnumOptions<Form> forms();
Expand Down
97 changes: 65 additions & 32 deletions src/math/function1D.cpp
View file
Expand Up @@ -67,8 +67,7 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
{Functions1D::Form::Gaussian,
{{"fwhm"},
{FunctionProperties::FourierTransform, FunctionProperties::Normalisation},
[](std::vector<double> params)
{
[](std::vector<double> params) {
params.push_back(params[0] / (2.0 * sqrt(2.0 * log(2.0))));
params.push_back(1.0 / params[1]);
return params;
Expand Down Expand Up @@ -105,8 +104,7 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
{Functions1D::Form::ScaledGaussian,
{{"A", "fwhm"},
{FunctionProperties::FourierTransform, FunctionProperties::Normalisation},
[](std::vector<double> params)
{
[](std::vector<double> params) {
params.push_back(params[1] / (2.0 * sqrt(2.0 * log(2.0))));
params.push_back(1.0 / params[2]);
return params;
Expand All @@ -116,17 +114,19 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
* f(x) = A exp ( - --------- )
* ( 2 * c * c )
*/
[](double x, double omega, const std::vector<double> &params)
{ return params[0] * exp(-(0.5 * x * x * params[3] * params[3])); },
[](double x, double omega, const std::vector<double> &params) {
return params[0] * exp(-(0.5 * x * x * params[3] * params[3]));
},
// First derivative (not defined)
{},
/*
* ( x * x * c * c )
* FT(x) = A exp ( - ------------- )
* ( 2 )
*/
[](double x, double omega, const std::vector<double> &params)
{ return params[0] * exp(-(0.5 * x * x * params[2] * params[2])); },
[](double x, double omega, const std::vector<double> &params) {
return params[0] * exp(-(0.5 * x * x * params[2] * params[2]));
},
/*
* 1
* Norm = --------------
Expand All @@ -144,8 +144,7 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
{Functions1D::Form::OmegaDependentGaussian,
{{"fwhm(x)"},
{FunctionProperties::FourierTransform, FunctionProperties::Normalisation},
[](std::vector<double> params)
{
[](std::vector<double> params) {
params.push_back(params[0] / (2.0 * sqrt(2.0 * log(2.0))));
params.push_back(1.0 / params[1]);
return params;
Expand All @@ -155,17 +154,19 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
* f(x) = exp ( - ---------------- )
* ( 2 * (c*omega)**2 )
*/
[](double x, double omega, const std::vector<double> &params)
{ return exp(-(x * x) / (2.0 * (params[1] * omega) * (params[1] * omega))); },
[](double x, double omega, const std::vector<double> &params) {
return exp(-(x * x) / (2.0 * (params[1] * omega) * (params[1] * omega)));
},
// First derivative (not defined)
{},
/*
* ( x*x * (c*omega)**2 )
* FT(x) = exp ( - ------------------ )
* ( 2 )
*/
[](double x, double omega, const std::vector<double> &params)
{ return exp(-(0.5 * x * x * (params[1] * omega) * (params[1] * omega))); },
[](double x, double omega, const std::vector<double> &params) {
return exp(-(0.5 * x * x * (params[1] * omega) * (params[1] * omega)));
},
/*
* 1
* Norm = ------------------
Expand All @@ -186,8 +187,7 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
{Functions1D::Form::GaussianC2,
{{"fwhm", "fwhm(x)"},
{FunctionProperties::FourierTransform},
[](std::vector<double> params)
{
[](std::vector<double> params) {
params.push_back(params[0] / (2.0 * sqrt(2.0 * log(2.0))));
params.push_back(params[1] / (2.0 * sqrt(2.0 * log(2.0))));
params.push_back(1.0 / params[2]);
Expand All @@ -199,24 +199,27 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
* f(x) = exp ( - ---------------------- )
* ( 2 * (c1 + c2*omega)**2 )
*/
[](double x, double omega, const std::vector<double> &params)
{ return exp(-(x * x) / (2.0 * (params[2] + params[3] * omega) * (params[2] + params[3] * omega))); },
[](double x, double omega, const std::vector<double> &params) {
return exp(-(x * x) / (2.0 * (params[2] + params[3] * omega) * (params[2] + params[3] * omega)));
},
// First derivative (not defined)
{},
/*
* ( x * x * (c1 + c2*omega)**2 )
* FT(x) = exp ( - -------------------------- )
* ( 2 )
*/
[](double x, double omega, const std::vector<double> &params)
{ return exp(-(0.5 * x * x * (params[2] + params[3] * omega) * (params[2] + params[3] * omega))); },
[](double x, double omega, const std::vector<double> &params) {
return exp(-(0.5 * x * x * (params[2] + params[3] * omega) * (params[2] + params[3] * omega)));
},
/*
* 1
* Norm = --------------------------
* (c1 + c2 omega) sqrt(2 pi)
*/
[](double omega, const std::vector<double> &params)
{ return 1.0 / ((params[2] + params[3] * omega) * sqrt(2.0 * M_PI)); }}},
[](double omega, const std::vector<double> &params) {
return 1.0 / ((params[2] + params[3] * omega) * sqrt(2.0 * M_PI));
}}},
/*
* Lennard-Jones 12-6 Potential
*
Expand All @@ -233,8 +236,7 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
* F(x) = 4 * epsilon * [ ( ----- ) - ( ----- ) ]
* [ ( x ) ( x ) ]
*/
[](double x, double omega, const std::vector<double> &params)
{
[](double x, double omega, const std::vector<double> &params) {
auto sigmar = params[1] / x;
auto sigmar6 = pow(sigmar, 6.0);
auto sigmar12 = sigmar6 * sigmar6;
Expand All @@ -245,25 +247,56 @@ static std::map<Functions1D::Form, Function1DDefinition> functions1D_ = {
* dYdX(x) = -48 * epsilon * [ ( --------- ) - 0.5 * ( -------- ) ]
* [ ( x**13 ) ( x**7 ) ]
*/
[](double x, double omega, const std::vector<double> &params)
{
[](double x, double omega, const std::vector<double> &params) {
auto sigmar = params[1] / x;
auto sigmar6 = pow(sigmar, 6.0);
return 48.0 * params[0] * sigmar6 * (-sigmar6 + 0.5) / x;
},
{},
{}}},
/*
* Buckingham Potential
*
* Parameters:
* INPUT 0 = A
* INPUT 1 = B
* INPUT 2 = C
*/
{Functions1D::Form::Buckingham,
{{"A", "B", "C"},
{FunctionProperties::FirstDerivative},
[](std::vector<double> params) { return params; },
/* C
* F(x)=A exp(-B * x) - -----
* x**6
*/
[](double x, double A, const std::vector<double> &params) {
auto B = exp(-params[1] * x);
auto C = params[2] * pow(x, 6.0);
return params[0] * B + C;
},
// dy/dx = -B * A exp(-B * x) + 6 * C * x**-7
[](double x, double A, const std::vector<double> &params) {
auto expo = exp(-params[1] * x);
auto C = 6 * params[2] * pow(x, -7.0);
return -params[1] * params[0] * expo + C;
},
{},
{}}}};

// Return enum option info for forms
EnumOptions<Functions1D::Form> Functions1D::forms()
{
return EnumOptions<Functions1D::Form>("Function1D",
{{Functions1D::Form::None, "None"},
{Functions1D::Form::Gaussian, "Gaussian", 1},
{Functions1D::Form::ScaledGaussian, "ScaledGaussian", 2},
{Functions1D::Form::OmegaDependentGaussian, "OmegaDependentGaussian", 1},
{Functions1D::Form::GaussianC2, "GaussianC2", 2},
{Functions1D::Form::LennardJones126, "LennardJones126", 2}});
{
{Functions1D::Form::None, "None"},
{Functions1D::Form::Gaussian, "Gaussian", 1},
{Functions1D::Form::ScaledGaussian, "ScaledGaussian", 2},
{Functions1D::Form::OmegaDependentGaussian, "OmegaDependentGaussian", 1},
{Functions1D::Form::GaussianC2, "GaussianC2", 2},
{Functions1D::Form::LennardJones126, "LennardJones126", 2},
{Functions1D::Form::Buckingham, "Buckingham", 3},
});
}

// Return parameters for specified form
Expand Down
3 changes: 2 additions & 1 deletion src/math/function1D.h
View file
Expand Up @@ -77,7 +77,8 @@ class Functions1D
ScaledGaussian,
OmegaDependentGaussian,
GaussianC2,
LennardJones126
LennardJones126,
Buckingham,
};
// Return enum options for form
static EnumOptions<Form> forms();
Expand Down

0 comments on commit 408996d

Please sign in to comment.