void qtrap(double a, double b, double& s)
{
    double eps = 0.000001;
    int jmax = 20;
    double olds = -1e+30;
    for (int j = 1; j<=jmax; j++)
	{
        trapzd(a, b, s, j);
        if (fabs(s - olds) < eps * fabs(olds))
		{
			return;
		}
        olds = s;
    }
    cout<<"too many steps."<<endl;
}