Language: C++

//
// Numerical Integration (Adaptive Gauss--Lobatto--Kronrod)
//
// Description:
//   Gauss--Lobatto formula is a numerical integrator
//   that is accurate for polynomials of degree <= 2n+1.
//   Adaptive Gauss--Lobatto recursively divides the domain
//   and computes integral by using G-L formula.
//
// Algorithm:
//   Above.
//
// Complexity:
//   O(#pieces) for a piecewise polynomials.
//   In general, it converges in O(1/n^6) for smooth functions.
//
// Verified: 
//   AOJ 2034
//
// References:
//   R. T. Trimbitas, M. G. Trimbitas (2008):
//   Gauss-Lobatto-Kronrod formulae and adaptive numerical integration.
//   in Proceedings of the 10th International Symposium on 
//   Symbolic and Numeric Algorithms for Scientific Computing. pp. 183--186.
//
#include <iostream>
#include <vector>
#include <limits>
#include <cstdio>
#include <algorithm>
#include <functional>

using namespace std;

#define fst first
#define snd second
#define all(c) ((c).begin()), ((c).end())

// adaptive gauss-lobatto
template <class F>
double integrate(F f, double lo, double hi, double eps = 1e-8) {
  const double th = eps / 1e-16; // (= eps / machine_epsilon)
  function<double (double,double,double,double)> rec = [&](double x0, double x6, double y0, double y6) {
    const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
    double x3 = (x0 + x6)/2, y3 = f(x3), h = (x6 - x0)/2;
    double x1 = x3-a*h, x2 = x3-b*h, x4 = x3+b*h, x5 = x3+a*h;
    double y1 = f(x1),  y2 = f(x2),  y4 = f(x4),  y5 = f(x5);
    double I1 = (y0+y6 + 5*(y2+y4)) * (h/6);
    double I2 = (77*(y0+y6) + 432*(y1+y5) + 625*(y2+y4) + 672*y3) * (h/1470);
    if (th + (I1 - I2) == th || x1 <= x0 || x5 >= x6) return I2;
    return (double)(rec(x0, x1, y0, y1) + rec(x1, x2, y1, y2) + rec(x2, x3, y2, y3) + 
                    rec(x3, x4, y3, y4) + rec(x4, x5, y4, y5) + rec(x5, x6, y5, y6));
  };
  return rec(lo, hi, f(lo), f(hi));
}

int main() {
  for (int n; scanf("%d", &n); ) {
    if (n == 0) break;
    double r; 
    scanf("%lf\n", &r);

    vector<pair<double,double>> p(n);
    for (int i = 0; i < n; ++i) 
      scanf("%lf %lf", &p[i].fst, &p[i].snd);
    function<double (double)> f = [&](double x) {
      if (x < p[0].fst) return 0.0;
      if (x >= p.back().fst) return 0.0;
      for (int i = 0; i+1 < n; ++i) 
        if (x < p[i+1].fst) 
          return p[i].snd + (p[i+1].snd-p[i].snd)/(p[i+1].fst-p[i].fst)*(x-p[i].fst);
    };
    function<double (double)> g = [&](double x) {
      return f(x) * f(x + r);
    };
    printf("%.12f\n", integrate(g, p[0].fst+min(0.0,r), p.back().fst+max(0.0,r)));
  }
}

Compiler output:
a.cc: In lambda function:
a.cc:44:22: error: there are no arguments to 'sqrt' that depend on a template parameter, so a declaration of 'sqrt' must be available [-fpermissive]
   44 |     const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
      |                      ^~~~
a.cc:44:22: note: (if you use '-fpermissive', G++ will accept your code, but allowing the use of an undeclared name is deprecated)
a.cc:44:47: error: there are no arguments to 'sqrt' that depend on a template parameter, so a declaration of 'sqrt' must be available [-fpermissive]
   44 |     const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
      |                                               ^~~~
a.cc: In instantiation of 'double integrate(F, double, double, double) [with F = std::function<double(double)>]':
a.cc:76:32:   required from here
   76 |     printf("%.12f\n", integrate(g, p[0].fst+min(0.0,r), p.back().fst+max(0.0,r)));
      |                       ~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
a.cc:44:26: error: 'sqrt' was not declared in this scope
   44 |     const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
      |                      ~~~~^~~~~~~~~
a.cc:44:51: error: 'sqrt' was not declared in this scope, and no declarations were found by argument-dependent lookup at the point of instantiation
   44 |     const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
      |                                               ~~~~^~~~~
a.cc:44:26: note: 'sqrt' declared here, later in the translation unit
   44 |     const double a = sqrt(2.0/3.0), b = 1.0 / sqrt(5.0);
      |                      ~~~~^~~~~~~~~
a.cc: In lambda function:
a.cc:72:5: warning: control reaches end of non-void function [-Wreturn-type]
   72 |     };
      |     ^