using System;
using static System.Console;
using static System.Double;
using static System.Math;
public static partial class matlib{

public static double open4
(Func<double,double> f,double a,double b,double acc=1e-3,double eps=1e-3,
double f2=NaN, double f3=NaN)
{
double h=b-a;
if(IsNaN(f2)){ f2=f(a+2*h/6); f3=f(a+4*h/6); } // first call, no points to reuse
double f1=f(a+h/6), f4=f(a+5*h/6);
double Q = (2*f1+f2+f3+2*f4)/6*(b-a);
double q = (  f1+f2+f3+  f4)/4*(b-a);
double err = Abs(Q-q);
double tol = acc+eps*Abs(Q);
if (err <= tol) return Q;
else return open4(f,a,(a+b)/2,acc/Sqrt(2),eps,f1,f2)
           +open4(f,(a+b)/2,b,acc/Sqrt(2),eps,f3,f4);
}

public static double open4cc
(Func<double,double> f,double a,double b,double acc=1e-3,double eps=1e-3)
{
Func<double,double> fcc = t => f((a+b)/2+(b-a)/2*Cos(t))*Sin(t)*(b-a)/2;
return open4(fcc,0,PI,acc,eps);
}

}//matlib