Where can I have an error in finding the area of ​​a complex figure?

Question

C++ / C#

Problem statement (be careful, large images)

I have written a solution but I am getting WA. Can you please tell me what could be wrong?

My solution code

#include <iostream>
#include <vector>
#include <algorithm>
#include <iomanip>
#include <fstream>
#include <cmath>

using namespace std;


class Sol {
public:

  vector<double> leftborders;
  vector<vector<double>> leftvalues;

  vector<double> rightborders;
  vector<vector<double>> rightvalues;

  double sum = 0;


  double calc(double a, double b, double c, double left, double right) {
    double r1 = a * pow(right, 3) / 3;
    double r2 = b * pow(right, 2) / 2;
    double r3 = c * right;


    double l1 = a * pow(left, 3) / 3;
    double l2 = b * pow(left, 2) / 2;
    double l3 = c * left;


    return r1 + r2 + r3 - (l1 + l2 + l3);

  }

  void doit() {

    double liters = 1;
    double riters = 1;

    while (liters <= leftvalues.size()) {

      double rside = min(leftborders[liters], rightborders[riters]);
      double lside = max(leftborders[liters-1], rightborders[riters-1]);

      double comp_a = leftvalues[liters-1][0] - rightvalues[riters-1][0];
      double comp_b = leftvalues[liters-1][1] - rightvalues[riters-1][1];
      double comp_c = leftvalues[liters-1][2] - rightvalues[riters-1][2];

      double desc = comp_b * comp_b - 4 * comp_a * comp_c;

      vector<double> stones;  

      if (desc == 0) {
        if (comp_a != 0) {
          double p1 = (-1 * comp_b) / (2 * comp_a);
          if (p1 > lside&& p1 < rside) stones.push_back(p1);
        }
      }
      else if (desc > 0) {
        if (comp_a == 0) {
          double p0 = (-1 * comp_c) / comp_b;
          if (p0 > lside&& p0 < rside) stones.push_back(p0);
        }
        else {
          double p1 = (-1 * comp_b - sqrt(desc)) / (2 * comp_a);
          if (p1 > lside&& p1 < rside) stones.push_back(p1);

          double p2 = (-1 * comp_b + sqrt(desc)) / (2 * comp_a);
          if (p2 > lside&& p2 < rside) stones.push_back(p2);
        }
      }
      stones.push_back(rside);

      for (int i = 0; i < stones.size(); ++i) {
        double one = calc(leftvalues[liters-1][0], leftvalues[liters-1][1], leftvalues[liters-1][2],lside,stones[i]);
        double two = calc(rightvalues[riters-1][0], rightvalues[riters-1][1], rightvalues[riters-1][2],lside,stones[i]);
        sum += abs(one-two);
        lside = stones[i];
      }

      if (leftborders[liters] <= rightborders[riters]) {
        liters++;
      }
      else {
        riters++;
      }
      

    }

  }

};


int main() {

  ifstream f("input.txt");
  Sol s;
  double leftorder;
  double rightorder;
  f >> leftorder;
  f >> rightorder;
  for (double i = 0; i <= leftorder; ++i) {
    double tmp;
    f >> tmp;
    s.leftborders.push_back(tmp);
  }

  for (double i = 0; i < leftorder; ++i) {
    double tmp1;
    double tmp2;
    double tmp3;
    f >> tmp1;
    f >> tmp2;
    f >> tmp3;
    vector<double> vals;
    vals.push_back(tmp1);
    vals.push_back(tmp2);
    vals.push_back(tmp3);
    s.leftvalues.push_back(vals);
  }

  for (double i = 0; i <= rightorder; ++i) {
    double tmp;
    f >> tmp;
    s.rightborders.push_back(tmp);
  }

  for (double i = 0; i < rightorder; ++i) {
    double tmp1;
    double tmp2;
    double tmp3;
    f >> tmp1;
    f >> tmp2;
    f >> tmp3;
    vector<double> vals;
    vals.push_back(tmp1);
    vals.push_back(tmp2);
    vals.push_back(tmp3);
    s.rightvalues.push_back(vals);
  }


  s.doit();
  cout << fixed << setprecision(10) << s.sum;

  return 0;
}

Answer 1

The discussion was in the comments, and the key reason is that when comp_a < 0, the roots are not sorted .
Side code issues.
1. It was possible to have an integral of a difference, but not a difference of integrals.
2. It's better not to use vector and pair instead of struct { a, b, c } - unless there is important automation.

Answer 2

Your code for merging two arrays of segments is not entirely clear. I would stupidly drive all intermediate points into one vector and sort it. Well, or made a standard merge of two sorted arrays. Further in the loop, I would count on the segment between two adjacent points (if the length of the segment is at least 1e-6), Curve parameters are easily searched - keep 2 counters, as you already have, and move them until the current segment for a particular function will close later than the start of the current segment. Maybe your code is equivalent to this, but I'm not sure.
Try to take some test and break the available pieces into several so that the two functions have completely different points, and so that they alternate irregularly (not the first, second, first, second ...). Try to have the points of the second curve be the last and vice versa.
Still, an obvious problem here comp_a != 0:. Floats can never be exactly compared. Epsilon only:

x < y  ====> x < y - eps
x <= y ====> x < y + eps 
x == y ====> fabs(x-y) < eps
x != y ====> fabs(x-y) > eps