#include <iostream>
#include <cmath>

using namespace std;

struct Point {
    double x;
    double y;

    Point() {x = y = 0;}
    Point(double _x, double _y) : x(_x), y(_y) {};
};

double df(Point x) {return 5 + x.x - x.y;}

void print(Point* x, int n) {
    cout << endl << "n = " << n << endl << endl;
    cout << "===== points =====" << endl;
    for (int i = 0; i < n + 1; i++) {
        cout << "[" << i << "]: " << "x = " << x[i].x << " y = " << x[i].y << endl;
        cout << endl;
    }
}

Point* rkSolve(Point s, double a, double b, double df(Point x), int n) {
    double k1, k2, k3, k4;
    double h = (b - a) / n;
    Point* res = new Point[n + 1];

    for (int i = 0; i < n + 1; i++) {
        res[i] = s;

        k1 = df(s);
        k2 = df(Point(s.x + 0.5 * h, s.y + 0.5 * h * k1));
        k3 = df(Point(s.x + 0.5 * h, s.y + 0.5 * h * k2));
        k4 = df(Point(s.x + h, s.y + h * k3));

        s.y = s.y + (h / 6) * (k1 + 2 * k2 + 2 * k3 + k4);
        s.x += h;
    }
    return res;
}

double checkEps(Point* s1, Point* s2, int n) {
    double eps = -1, curEps;

    for (int i = 1; i < n + 1; i++) {
        curEps = fabs(s1[i].y - s2[i * 2].y) / 15;
        eps = max(eps, curEps);
    }
    return eps;
}

void solve(Point s, double a, double b, double df(Point x), double eps) {
    int n = 10;
    Point *s1, *s2;

    do {
        n *= 2;

        s1 = rkSolve(s, a, b, df, n);
        s2 = rkSolve(s, a, b, df, n * 2);

    } while (checkEps(s1, s2, n) >= eps);

    print(s1, n);
}

int main() {
    Point s;
    double a, b, eps;

    cout << "Enter start point: ";
    cin >> s.x >> s.y;

    cout << "Enter a & b: ";
    cin >> a >> b;

    cout << "Enter eps: ";
    cin >> eps;

    solve(s, a, b, df, eps);

    return 0;
}