Find the area between two functions

I solve the following problem.

There are two functions f(x) and g(x). Each of them is divided into parts. The f(x) function consists of 'n' parts, and g(x) consists of 'm' parts. Each part is represented as a Second-Degree Polynomial Function (ax^2 + bx + c). We need to find an area between f(x) and g(x) with precision 10^(-6).

Input data:

n, m (1 <= n, m <= 10^5);

boundary points of f(x); //n+1 points

polynomials of f(x); //n strings (a, b, c for ax^2 + bx + c)

boundary points of g(x); //m+1 points

polynomials of g(x); //m strings (a, b, c for ax^2 + bx + c)

The first and last points of the functions are equal.

Sample input:

Output:

1.3333333333

Here is my code:

#include <bits/stdc++.h>
using namespace std;


long double integrate(int f[3], int g[3], long double left, long double right) {
    long double a = f[0] - g[0], b = f[1] - g[1], c = f[2] - g[2];
    long double up = (a * right * right * right) / 3 + (b * right * right) / 2 + c * right;
    long double down = (a * left * left * left) / 3 + (b * left * left) / 2 + c * left;
    return (up > down) ? up - down : down - up;
}

void solve(int f[3], int g[3], int left, int right, long double &sum) {

    long double a = f[0] - g[0], b = f[1] - g[1], c = f[2] - g[2];
    if (a == 0) {
        if (b != 0) {
            long double x = -c/b;
            if (x > left && x < right) {
                sum += integrate(f, g, left, x);
                sum += integrate(f, g, x, right);
            }
        } else {
            sum += integrate(f, g, left, right);
        }
        return;
    }

    long double discriminant = b * b - 4 * a * c;
    if (discriminant < 0) { sum += integrate(f, g, left, right); }
    else {
        long double q = b >= 0 ? (-b - sqrt(discriminant))/2 : (-b + sqrt(discriminant))/2;
        long double x1 = q / a, x2 = c / q;
        if (discriminant == 0.0) {
            if (x1 > left && x1 < right) {
                sum += integrate(f, g, left, x1);
                sum += integrate(f, g, x1, right);
            } else {
                sum += integrate(f, g, left, right);
            }
        } else {
            long double first = min(x1, x2), second = max(x1, x2);
            if (first > left && first < right) {
                sum += integrate(f, g, left, first);
                if (second > left && second < right) {
                    sum += integrate(f, g, first, second);
                    sum += integrate(f, g, second, right);
                } else {
                    sum += integrate(f, g, first, right);
                }
            } else if (second > left && second < right) {
                sum += integrate(f, g, left, second);
                sum += integrate(f, g, second, right);
            } else {
                sum += integrate(f, g, left, right);
            }
        }
    }
    return;
}


int main() {
//freopen("input.txt", "r", stdin);
//freopen("output.txt", "w", stdout);

    int n, m;
    cin >> n >> m;

    int f[n+1];
    int ffunc[n][3];
    int g[m+1];
    int gfunc[m][3];
    set <int> points;

    for (int i = 0; i < n+1; i++) {
        cin >> f[i];
        points.insert(f[i]);
    }

    for (int i = 0; i < n; i++) {
        for (int k = 0; k < 3; k++) {
            cin >> ffunc[i][k];
        }
    }

    for (int i = 0; i < m+1; i++) {
        cin >> g[i];
        points.insert(g[i]);
    }

    for (int i = 0; i < n; i++) {
        for (int k = 0; k < 3; k++) {
            cin >> gfunc[i][k];
        }
    }

    int fit = 0, git = 0;
    long double sum = 0.0;
    auto it1 = points.begin();
    auto it2 = points.begin();
    it2++;
    while (it2 != points.end()) {

        solve(ffunc[fit], gfunc[git], *it1, *it2, sum);

        if (f[fit+1] == *it2) {
            fit++;
        }
        if (g[git+1] == *it2) {
            git++;
        }
        it1++;
        it2++;
    }
    cout.precision(27);
    cout << sum;
    return 0;
}

The program doesn't pass some tests. It gets wrong answer.

What could be the problem?

Solution

In solve, there is at least a corner case not accounted for:

if (a == 0) {
    if (b != 0) {
        long double x = -c/b;
        if (x > left && x < right) {
            sum += integrate(f, g, left, x);
            sum += integrate(f, g, x, right);
        }  
        // else {
        //     sum += integrate(f, g, left, right);
        // }
    } else {
        sum += integrate(f, g, left, right);
    }
    return;
}

In main, the local variables f, ffunc, g and gfunc are all variable length arrays (a non-standard extension provided by some compilers) and should rather all be standard containers like std::vector.

The posted code also uses a std::set to store all the points and to ensure their order. This may not be the most efficient way to accomplish that result, if the boundary points are already ordered in their own ranges.

Here, a possible different implementation is tested.