There are advanced and multiple variants of Gaussian quadrature formula which give good results already.

It is necessary to specify the word "optimal".

It's a naive but an interesting question. Naive in the sense, as Vladislav pointed out, one need to be more specific: for example:

• is it for multivariate function (integrand)
• regularity of the integrand, for example, one need to see the singularity. etc etc etc....

A straightforward answer to my knowledge is the answer is "problem-specific". Nevertheless, Gauss quadrature based on orthogonal polynomials... seems to be reasonable (don't hang me if it won't work for you!). You may have a look at Nick Trefethen's book, "Spectral methods in Matlab".

P.S.: An elaborate discussion on the similar question you can find at this link: https://www.researchgate.net/post/What_is_the_best_Numerical_integration_method_and_why

https://www.researchgate.net/post/What_is_the_best_Numerical_integration_method_and_why

If you want to second order accuracy, the formula of average rectangles:

SUMMA{f((xi+1+xi)/2) * (xi+1-xi)}

more precisely than the formula of trapezium:

SUMMA{ (f(xi+1)+f(xi) )/2 * (xi+1-xi)}

This is especially visible when the integrand expression f(xi) obtained numerically.

 You can get the optimal quadrature integration, by using determine coefficients method, which we get the accuracy from truncation errors.  

I think the exhaustive discussion is under the link above. If You do not like to study all aspects of problem and theory (I think it is not good idea) the simplest way to use simplest method and then to use Simpson rule. I/e/

You can calculate the integral for some set of discretization steps and then by some estimation of error You can obtain corrected solution.

For example Y1/2=Y+a*(1/2)^n

Y1/4=Y+a*(1/4)^n

Y1/8=Y+a*(1/8)^n

If you calculate three numerical results You can determine a,n and Y.

For example, it works well for integral from sqrt(x) on interval (0,1).

Thank you very much for the answers to my question!

I've searched for some DEFINITION for the term "optimal quadrature" and an ALGORITHM to construct such a formula for specified class of functions. Seems that this term is not properly defined.

Vladislav Galkin, Tapan Kumar Hota, you are absolutely right that my question is naive in that way  and It is necessary to specify the problem. So, I give some explanation below.

My question arised while I was reading the article by S.M.Sadatrasoul and M.Ezzati "Numerical solution of two-dimensional nonlinear Hammerstein fuzzy integral equation based on optimal fuzzy quadrature formula", Journal of Computational and Applied Mathematics 292 (2016) 430-446.

I want to apply their formulas for another kind of integral equations. Or to construct something similar. But I wonder what they mean using the term "optimal formula". I investigate one-dimensional case firtst. So, you can find in Theorem 19 formula (4) - the authors  say that it is the OPTIMAL quadrature formula among all formulas (3). For formula (3) they say that it has the minimal error among all quadrature formulas that use given knots. If the knots are given (fixed) we cannot use Gaussian quadrature formula as Vikash Pandey and Tapan Kumar Hota suggest, because in Gaussian formulas usually we use the roots of some special orthogonal polynomials (we define the knots). Of course, I can try different Gaussian quadrature for different kind of functions (polynomial, exponential, trigonometric). But it seem that it is not the way S.M.Sadatrasoul and M.Ezzati did it...

Vladislav Galkin, you show the average rectangles formulas - it is the same formula as (4) in cited above paper... I've already made a comparison between trapezoid and that formula and it is seen that their "optimal" formula gives two times smaller error... Of course, it is a known fact from the textbooks.

Important feature that is numerically verified. If xi irregular and values functions  f(xi) contain noise from the calculations, the average rectangles formulas reduces the noise better than trapezoid. This is important when after calculating fluid parameters, flow consumption must be found.

"Best" depends on the problem. GQ is the best all-round method. There is a zip file on my web site that contains every quadrature method you could possibly want, ready to cut-and-paste into your project.

A best way for solving multivariate/high dimensional integrations can by solved by the method of Laplace approximation or saddle-point approximation, however these methods are complicated. For simple concepts please see : https://www.hindawi.com/journals/aaa/2013/426916/ and http://ads.harvard.edu/books/1990fnmd.book/chapt5.pdf