bonds.h

#pragma once
#include "particle.h"

namespace Faunus::pairpotential {

struct BondData
{
    enum class Variant
    {
        HARMONIC = 0,
        FENE,
        FENEWCA,
        HARMONIC_TORSION,
        GROMOS_TORSION,
        PERIODIC_DIHEDRAL,
        HARMONIC_DIHEDRAL,
        INVALID
    };
    std::vector<int> indices; 

    std::function<double(Geometry::DistanceFunction)> energyFunc = nullptr;

    using IndexAndForce =
        std::pair<int, Point>; 

    std::function<std::vector<IndexAndForce>(Geometry::DistanceFunction)> forceFunc = nullptr;

    virtual void from_json(const json&) = 0;
    virtual void to_json(json&) const = 0;
    [[nodiscard]] virtual int numindex() const = 0; 
    [[nodiscard]] virtual Variant type() const = 0; 
    [[nodiscard]] virtual std::shared_ptr<BondData>
    clone() const = 0; 
    virtual void setEnergyFunction(
        const ParticleVector& particles) = 0;     
    [[nodiscard]] bool hasEnergyFunction() const; 
    [[nodiscard]] bool hasForceFunction() const;  
    void shiftIndices(int offset);                
    BondData() = default;
    explicit BondData(const std::vector<int>& indices);
    virtual ~BondData() = default;
};

NLOHMANN_JSON_SERIALIZE_ENUM(BondData::Variant,
                             {{BondData::Variant::INVALID, nullptr},
                              {BondData::Variant::HARMONIC, "harmonic"},
                              {BondData::Variant::FENE, "fene"},
                              {BondData::Variant::FENEWCA, "fene+wca"},
                              {BondData::Variant::HARMONIC_TORSION, "harmonic_torsion"},
                              {BondData::Variant::GROMOS_TORSION, "gromos_torsion"},
                              {BondData::Variant::PERIODIC_DIHEDRAL, "periodic_dihedral"},
                              {BondData::Variant::HARMONIC_DIHEDRAL, "harmonic_dihedral"}})

struct StretchData : public BondData
{
    [[nodiscard]] int numindex() const override;
    StretchData() = default;
    explicit StretchData(const std::vector<int>& indices);
};

struct TorsionData : public BondData
{
    [[nodiscard]] int numindex() const override;
    TorsionData() = default;
    explicit TorsionData(const std::vector<int>& indices);
};

struct HarmonicBond : public StretchData
{
    double half_force_constant = 0.0;
    double equilibrium_distance = 0.0;
    [[nodiscard]] Variant type() const override;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    void
    setEnergyFunction(const ParticleVector& particles) override; 
    HarmonicBond() = default;
    HarmonicBond(double k, double req, const std::vector<int>& indices);
};

struct FENEBond : public StretchData
{
    double half_force_constant = 0.0;
    double max_squared_distance = 0.0;
    [[nodiscard]] Variant type() const override;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    void setEnergyFunction(const ParticleVector& particles) override;
    FENEBond() = default;
    FENEBond(double k, double rmax, const std::vector<int>& indices);
};

struct FENEWCABond : public StretchData
{
    double half_force_constant = 0.0;
    double max_distance_squared = 0.0;
    double epsilon = 0.0;
    double sigma_squared = 0.0;
    std::array<double, 4> k = {{0.0, 0.0, 0.0, 0.0}};
    [[nodiscard]] Variant type() const override;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    void setEnergyFunction(const ParticleVector& calculateDistance) override;
    FENEWCABond() = default;
    FENEWCABond(double k, double rmax, double epsilon, double sigma,
                const std::vector<int>& indices);
};

struct HarmonicTorsion : public TorsionData
{
    double half_force_constant = 0.0;
    double equilibrium_angle = 0.0;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    [[nodiscard]] Variant type() const override;
    void setEnergyFunction(const ParticleVector& particles) override;
    HarmonicTorsion() = default;
    HarmonicTorsion(double k, double aeq, const std::vector<int>& indices);
};

struct GromosTorsion : public TorsionData
{
    double half_force_constant = 0.0;
    double cosine_equilibrium_angle = 0.0;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    [[nodiscard]] Variant type() const override;
    void setEnergyFunction(const ParticleVector& calculateDistance) override;
    GromosTorsion() = default;
    GromosTorsion(double k, double cos_aeq, const std::vector<int>& indices);
};

struct PeriodicDihedral : public BondData
{
    double force_constant = 0.0;
    double phase_angle = 0.0;
    double periodicity = 1.0;
    [[nodiscard]] int numindex() const override;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    [[nodiscard]] Variant type() const override;
    void setEnergyFunction(const ParticleVector& particles) override;
    PeriodicDihedral() = default;
    PeriodicDihedral(double k, double phi, double n, const std::vector<int>& indices);
};

struct HarmonicDihedral : public BondData
{
    double half_force_constant = 0.0;
    double equilibrium_dihedral = 0.0;
    [[nodiscard]] int numindex() const override;
    [[nodiscard]] std::shared_ptr<BondData> clone() const override;
    void from_json(const json& j) override;
    void to_json(json& j) const override;
    [[nodiscard]] Variant type() const override;
    void setEnergyFunction(const ParticleVector& particles) override;
    HarmonicDihedral() = default;
    HarmonicDihedral(double k, double deq, const std::vector<int>& indices);
};

void from_json(const json& j, std::shared_ptr<BondData>& bond);
void to_json(json& j, const std::shared_ptr<const BondData>& bond);
void to_json(json& j, const BondData& bond);

} // namespace Faunus::pairpotential