CoreFunction.h

Go to the documentation of this file.

// n2p2 - A neural network potential package
// Copyright (C) 2018 Andreas Singraber (University of Vienna)
// Copyright (C) 2020 Martin P. Bircher
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see <https://www.gnu.org/licenses/>.
 
#ifndef COREFUNCTION_H
#define COREFUNCTION_H
 
#include <cmath>
#include <string>
#include <vector>
 
namespace nnp
{
 
class CoreFunction
{
public:
    enum class Type
    {
        COS,
        POLY1,
        POLY2,
        POLY3,
        POLY4,
        EXP
    };
 
    CoreFunction();
    void                     setType(Type const type);
    void                     setType(std::string const typeString);
    Type                     getType() const;
#ifndef N2P2_NO_ASYM_POLY
    void                     setAsymmetric(bool asymmetric);
    bool                     getAsymmetric() const;
#endif
    double                   f(double x) const;
    double                   df(double x) const;
    void                     fdf(double x, double& fx, double& dfx) const;
    std::vector<std::string> info() const;
 
private:
    static double const PI;
    static double const PI_2;
    static double const E;
    Type                type;
#ifndef N2P2_NO_ASYM_POLY
    bool                asymmetric;
#endif
    double     (CoreFunction::*fPtr)(double x) const;
    double     (CoreFunction::*dfPtr)(double x) const;
    void       (CoreFunction::*fdfPtr)(double  x,
                                       double& fx,
                                       double& dfx) const;
 
    // Individual functions.
    double   fCOS(double x) const;
    double  dfCOS(double x) const;
    void   fdfCOS(double x, double& fx, double& dfx) const;
 
    double   fPOLY1(double x) const;
    double  dfPOLY1(double x) const;
    void   fdfPOLY1(double x, double& fx, double& dfx) const;
 
    double   fPOLY2(double x) const;
    double  dfPOLY2(double x) const;
    void   fdfPOLY2(double x, double& fx, double& dfx) const;
 
    double   fPOLY3(double x) const;
    double  dfPOLY3(double x) const;
    void   fdfPOLY3(double x, double& fx, double& dfx) const;
 
    double   fPOLY4(double x) const;
    double  dfPOLY4(double x) const;
    void   fdfPOLY4(double x, double& fx, double& dfx) const;
 
    double   fEXP(double x) const;
    double  dfEXP(double x) const;
    void   fdfEXP(double x, double& fx, double& dfx) const;
};
 
// Inlined function definitions //
 
inline CoreFunction::Type CoreFunction::getType() const { return type; }
 
#ifndef N2P2_NO_ASYM_POLY
inline void CoreFunction::setAsymmetric(bool asymmetric)
{
    this->asymmetric = asymmetric;
 
    return;
}
 
inline bool CoreFunction::getAsymmetric() const { return asymmetric; }
#endif
 
inline double CoreFunction::f(double x) const
{
#ifndef N2P2_NO_ASYM_POLY
    if (asymmetric) x = (2.0 - x) * x;
#endif
    return (this->*fPtr)(x);
}
 
inline double CoreFunction::df(double x) const
{
#ifndef N2P2_NO_ASYM_POLY
    if (asymmetric) return (2.0 - 2.0 * x) * (this->*dfPtr)((2.0 - x) * x);
    else return (this->*dfPtr)(x);
#else
    return (this->*dfPtr)(x);
#endif
}
 
inline void CoreFunction::fdf(double x, double& fx, double& dfx) const
{
#ifndef N2P2_NO_ASYM_POLY
    if (asymmetric)
    {
        (this->*fdfPtr)((2.0 - x) * x, fx, dfx);
        dfx *= 2.0 - 2.0 * x;
    }
    else (this->*fdfPtr)(x, fx, dfx);
#else
    (this->*fdfPtr)(x, fx, dfx);
#endif
    return;
}
 
inline double CoreFunction::fCOS(double x) const
{
    return 0.5 * (cos(PI * x) + 1.0);
}
 
inline double CoreFunction::dfCOS(double x) const
{
    return -PI_2 * sin(PI * x);
}
 
inline void CoreFunction::fdfCOS(double x, double& fx, double& dfx) const
{
    double const temp = cos(PI * x);
    fx = 0.5 * (temp + 1.0);
    dfx = -0.5 * PI * sqrt(1.0 - temp * temp);
    return;
}
 
inline double CoreFunction::fPOLY1(double x) const
{
    return (2.0 * x - 3.0) * x * x + 1.0;
}
 
inline double CoreFunction::dfPOLY1(double x) const
{
    return x * (6.0 * x - 6.0);
}
 
inline void CoreFunction::fdfPOLY1(double x, double& fx, double& dfx) const
{
    fx = (2.0 * x - 3.0) * x * x + 1.0;
    dfx = x * (6.0 * x - 6.0);
    return;
}
 
inline double CoreFunction::fPOLY2(double x) const
{
    return ((15.0 - 6.0 * x) * x - 10.0) * x * x * x + 1.0;
}
 
inline double CoreFunction::dfPOLY2(double x) const
{
    return x * x * ((60.0 - 30.0 * x) * x - 30.0);
}
 
inline void CoreFunction::fdfPOLY2(double x, double& fx, double& dfx) const
{
    double const x2 = x * x;
    fx = ((15.0 - 6.0 * x) * x - 10.0) * x * x2 + 1.0;
    dfx = x2 * ((60.0 - 30.0 * x) * x - 30.0);
    return;
}
 
inline double CoreFunction::fPOLY3(double x) const
{
    double const x2 = x * x;
    return (x * (x * (20.0 * x - 70.0) + 84.0) - 35.0) * x2 * x2 + 1.0;
}
 
inline double CoreFunction::dfPOLY3(double x) const
{
    return x * x * x * (x * (x * (140.0 * x - 420.0) + 420.0) - 140.0);
}
 
inline void CoreFunction::fdfPOLY3(double x, double& fx, double& dfx) const
{
    double const x2 = x * x;
    fx = (x * (x * (20.0 * x - 70.0) + 84.0) - 35.0) * x2 * x2 + 1.0;
    dfx = x2 * x * (x * (x * (140.0 * x - 420.0) + 420.0) - 140.0);
    return;
}
 
inline double CoreFunction::fPOLY4(double x) const
{
    double const x2 = x * x;
    return (x * (x * ((315.0 - 70.0 * x) * x - 540.0) + 420.0) - 126.0) *
           x2 * x2 * x + 1.0;
}
 
inline double CoreFunction::dfPOLY4(double x) const
{
    double const x2 = x * x;
    return x2 * x2 *
           (x * (x * ((2520.0 - 630.0 * x) * x - 3780.0) + 2520.0) - 630.0);
}
 
inline void CoreFunction::fdfPOLY4(double x, double& fx, double& dfx) const
{
    double x4 = x * x;
    x4 *= x4;
    fx = (x * (x * ((315.0 - 70.0 * x) * x - 540.0) + 420.0) - 126.0) *
         x * x4 + 1.0;
    dfx = x4 *
          (x * (x * ((2520.0 - 630.0 * x) * x - 3780.0) + 2520.0) - 630.0);
    return;
}
 
inline double CoreFunction::fEXP(double x) const
{
    return E * exp(1.0 / (x * x - 1.0));
}
 
inline double CoreFunction::dfEXP(double x) const
{
    double const temp = 1.0 / (x * x - 1.0);
    return -2.0 * E * x * temp * temp * exp(temp);
}
 
inline void CoreFunction::fdfEXP(double x, double& fc, double& dfc) const
{
    double const temp = 1.0 / (x * x - 1.0);
    double const temp2 = exp(temp);
    fc = E * temp2;
    dfc = -2.0 * E * x * temp * temp * temp2;
    return;
}
 
}
 
#endif