#include <iostream>
#include <vector>
#include <random>
#include <algorithm>
#include <cmath>
#include <limits>

inline int32_t GenerateRandomInt(int32_t lowerBound, int32_t upperBound)
{
    static std::random_device rd;
    static std::mt19937 prng(rd());
    std::uniform_int_distribution distrib(lowerBound, upperBound);
    return distrib(prng);
}

struct Point
{
    double x;
    double y;
};

template <auto Foo>
class ImmuneSystem
{
    class Lymphocyte
    {
    public:
        Lymphocyte() = default;

        void Mutate()
        {
            GenerateRandomInt(0, 1) ? MutateA() : MutateB();
        }

        std::pair<double, double> GetCoeffs()
        {
            return { m_a, m_b };
        }

        void RecalculateFitnessMetric(const std::vector<Point>& points)
        {
            double value = 0;

            for (const auto& [x, y] : points)
            {
                value += std::pow(y - Foo(m_a, m_b, x), 2);
            }

            m_fitness = std::sqrt(value);
        }

        double FitnessMetric() const
        {
            return m_fitness;
        }

    private:

        double Mutate(double coeff) const
        {
            return std::min(std::max(MinCoeffValue, coeff + GenerateRandomInt(-MutateBound, MutateBound)), MaxCoeffValue);
        }

        void MutateA()
        {
            m_a = Mutate(m_a);
        }

        void MutateB()
        {
            m_b = Mutate(m_b);
        }

    private:
        static constexpr double MutateBound = 5;
        static constexpr double MaxCoeffValue = 100;
        static constexpr double MinCoeffValue = -100;

        double m_a = 0.0;
        double m_b = 0.0;
        double m_fitness = 0;
    };

public:
    explicit ImmuneSystem(const std::vector<Point>& points, size_t iterationsCount, size_t lymphocytesCount)
    {
        std::vector<Lymphocyte> lymphocytes(lymphocytesCount);

        const auto comparator = [](const Lymphocyte& first, const Lymphocyte& second)
        {
            return first.FitnessMetric() < second.FitnessMetric();
        };

        for (auto& lymphocyte : lymphocytes)
            lymphocyte.RecalculateFitnessMetric(points);

        for (size_t i = 0; i < iterationsCount; ++i)
        {
            std::cout << "Iter: " << i << std::endl;
            std::sort(lymphocytes.begin(), lymphocytes.end(), comparator);

            for (size_t j = 0; j < lymphocytes.size() / 2; ++j)
            {
                lymphocytes[j].Mutate();
            }

            for (auto& lymphocyte : lymphocytes)
                lymphocyte.RecalculateFitnessMetric(points);
        }

        m_bestLymphocyte = lymphocytes.front();
    }

    std::pair<double, double> GetCoefficients()
    {
        return m_bestLymphocyte.GetCoeffs();
    }

    Lymphocyte m_bestLymphocyte;
};

template <auto Foo>
std::vector<Point> GeneratePoints(double a, double b, size_t count)
{
    std::vector<Point> points;
    points.reserve(count);
    for (size_t i = 0; i < count; ++i)
    {
        double x = static_cast<double>(GenerateRandomInt(std::numeric_limits<int16_t>::min(), std::numeric_limits<int16_t>::max()));
        points.push_back({ x, Foo(a, b, x) });
    }

    return points;
}

int main()
{
    constexpr auto func = [](double a, double b, double x)
    {
        return a + b / x;
    };

    static constexpr size_t IterationsCount = 500;
    static constexpr size_t LymphocytesCount = 100;
    static constexpr size_t PointsCount = 1000;

    static constexpr double A = -12.3;
    static constexpr double B = 45.3;

    ImmuneSystem<func> is(GeneratePoints<func>(A, B, PointsCount), IterationsCount, LymphocytesCount);
    const auto&& [a, b] = is.GetCoefficients();

    std::cout << "Coef A: " << a << ", Coef B: " << b << std::endl;

    return 0;
}