{"id":1576,"date":"2016-05-23T12:00:32","date_gmt":"2016-05-23T17:00:32","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=1576"},"modified":"2016-12-05T14:07:18","modified_gmt":"2016-12-05T19:07:18","slug":"modernization-of-numerical-integration-from-quad-to-integral","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2016\/05\/23\/modernization-of-numerical-integration-from-quad-to-integral\/","title":{"rendered":"Modernization of Numerical Integration, From Quad to Integral"},"content":{"rendered":"<div class=\"content\"><!--introduction-->The MATLAB functions for the numerical evaluation of integrals has evolved from <tt>quad<\/tt>, through <tt>quadl<\/tt> and <tt>quadgk<\/tt>, to today's <tt>integral<\/tt>.\r\n\r\n<!--\/introduction-->\r\n<h3>Contents<\/h3>\r\n<div>\r\n<ul>\r\n \t<li><a href=\"#b5f8b6f9-2641-409e-a12f-033dcab17486\">Quadrature<\/a><\/li>\r\n \t<li><a href=\"#2fc72a13-ca0d-4201-8dd1-e13bc45b12a6\">Adaptive Simpson's method and <tt>quad<\/tt><\/a><\/li>\r\n \t<li><a href=\"#3213406a-41cb-4e21-b1f6-46b5ace7384d\"><tt>Humps<\/tt> and <tt>quadgui<\/tt><\/a><\/li>\r\n \t<li><a href=\"#a88a5eb2-8cee-4496-a655-2416c5dae7f7\">Lobatto, Kronrod and <tt>quadl<\/tt><\/a><\/li>\r\n \t<li><a href=\"#07a23423-d0df-4c19-a51e-3c163ad0b85f\">Gauss, Kronrod, and <tt>quadgk<\/tt><\/a><\/li>\r\n \t<li><a href=\"#8f394fe9-1306-4888-ba9c-a2d792dad285\">Comparison on <tt>humps<\/tt>.<\/a><\/li>\r\n \t<li><a href=\"#b20b2fd4-8c24-48a3-bfb8-804976454356\"><tt>quadgk<\/tt> becomes <tt>integral<\/tt><\/a><\/li>\r\n \t<li><a href=\"#dc1fbbc4-e550-4971-9509-1fb0078e1939\">Reference<\/a><\/li>\r\n<\/ul>\r\n<\/div>\r\n<h4>Quadrature<a name=\"b5f8b6f9-2641-409e-a12f-033dcab17486\"><\/a><\/h4>\r\nSeveral years ago I regarded myself as an expert on numerical methods for approximating integrals. But then I stopped paying attention for a while. So this is an opportunity to get reacquainted. The task involves a real-valued function $f(x)$ of a real variable $x$, defined on an interval $a \\le x \\le b$. We seek to compute the value of the definite integral,\r\n\r\n$$ \\int_a^b \\ f(x) \\ dx $$\r\n\r\nI wrote a chapter on the subject in my textbook <i>Numerical Computing with MATLAB<\/i>. Here is <a href=\"https:\/\/www.mathworks.com\/content\/dam\/mathworks\/mathworks-dot-com\/moler\/quad.pdf\">a link to that chapter<\/a>. The chapter is titled \"Quadrature\". One of the first things I have to do in that chapter is explain that <i>quadrature<\/i> is an old-fashioned term derived from the notion that the value of a definite integral can be approximated by plotting the function on graph paper and then counting the number of little squares that fall underneath the graph. Here is a coarse picture.\r\n\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/quads.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n<h4>Adaptive Simpson's method and <tt>quad<\/tt><a name=\"2fc72a13-ca0d-4201-8dd1-e13bc45b12a6\"><\/a><\/h4>\r\nFor most of the time since the beginning of MATLAB, up until when I wrote the chapter on quadrature twelve years ago, the MATLAB function for evaluating definite integrals was one named <tt>quad<\/tt>. It was based on adaptive Simpson's method, pioneered by my grad school buddy Bill McKeeman to show off recursion in Algol.\r\n\r\nThese two figures illustrate the basis for adaptive Simpson's method.\r\n\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/simpsons.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\nThe algorithm works recursively on subintervals of the original interval [a,b]. If the values obtained from the 3-point Simpson's rule and the 5-point composite Simpson's rule agree to within a specified tolerance (which they clearly do not for these figures), then an extrapolation taken from the two values is accepted as the value for the integral over the subinterval. If not, the algorithm is recursively applied to the two halves of the subinterval.\r\n<h4><tt>Humps<\/tt> and <tt>quadgui<\/tt><a name=\"3213406a-41cb-4e21-b1f6-46b5ace7384d\"><\/a><\/h4>\r\nI have used a function named <tt>humps<\/tt> to demonstrate numerical integration for a long time. It has a strong peak at x = 0.3 and a milder peak at x = 0.9. It has found its way into the MATLAB <tt>demos<\/tt> directory, so you have it on your own machine. A graph appears below.\r\n\r\n$$ h(x) = \\frac{1}{(x-0.3)^2+.01} + \\frac{1}{(x-0.9)^2+.04} - 6 $$\r\n\r\nAn interactive graphical app named <tt>quadgui<\/tt> is included with the programs for <i>Numerical Computing with MATLAB<\/i>, available from <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/37976-numerical-computing-with-matlab\">MATLAB Central File Exchange<\/a>. It shows <tt>quad<\/tt> in action. Here is a snapshot of <tt>quad<\/tt> part way through computing the integral of <tt>humps<\/tt> over the interval [0,1]. Integration over the interval to the left of 1\/4 is finished. The integration over the subinterval [1\/4,3\/8] has just been successfully completed. A few samples to the right of 3\/8 were computed during the initialization. The integrand has been evaluated 19 times.\r\n\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/step19.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\nHere is the final result. The default tolerance for <tt>quad<\/tt> is $10^{-6}$ and the value of the integral obtained to six significant digits is 29.8583. The integrand has been evaluated 93 times. Most of those evaluations took place near the peaks where the <tt>humps<\/tt> varies rapidly. This indicates that the adaptivity works as expected.\r\n\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/humps.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n\r\n<tt>Humps<\/tt> is a rational function. It is possible to use the Symbolic Toolbox to obtain the indefinite and then the definite integral analytically.\r\n<pre class=\"codeinput\">   syms <span class=\"string\">x<\/span>\r\n   h = 1\/((x-3\/10)^2+1\/100) + 1\/((x-9\/10)^2+4\/100) - 6\r\n   I = int(h,x)\r\n   D = subs(I,1)-subs(I,0)\r\n   format <span class=\"string\">long<\/span>\r\n   Q = double(D)\r\n<\/pre>\r\n<pre class=\"codeoutput\">h =\r\n1\/((x - 9\/10)^2 + 1\/25) + 1\/((x - 3\/10)^2 + 1\/100) - 6\r\nI =\r\n10*atan(10*x - 3) - 6*x + 5*atan(5*x - 9\/2)\r\nD =\r\n5*atan(1\/2) + 10*atan(3) + 10*atan(7) + 5*atan(9\/2) - 6\r\nQ =\r\n  29.858325395498674\r\n<\/pre>\r\nThis confirms the six figure value shown by <tt>quadqui<\/tt>\r\n<h4>Lobatto, Kronrod and <tt>quadl<\/tt><a name=\"a88a5eb2-8cee-4496-a655-2416c5dae7f7\"><\/a><\/h4>\r\nAbout the time in 2004 that I was writing the <i>Quadrature<\/i> chapter for <i>NCM<\/i>, the <tt>quadl<\/tt> function showed up. This uses the <i>Lobatto<\/i> rule.\r\n\r\n$$ \\int_{-1}^1 \\ f(x) \\ dx \\ \\approx \\ w_1 f(-1) + w_2 f(-x_1)\r\n+ w_2 f(x_1) + w_1 f(1) $$\r\n\r\nwith $x_1 = 1\/\\sqrt{5}, \\ w_1 = 1\/6, \\ w_2 = 5\/6$.\r\n\r\n<tt>Quadl<\/tt> also employs the associated <i>Kronrod<\/i> rule, which provides 13 additional irrational abscissa and weights. The Lobatto-Kronrod rule is much higher order and consequently more accurate than Simpson's rule. As a result <tt>quadl<\/tt> is usually faster and more accurate than <tt>quad<\/tt>. But I just devoted a couple of exercises in the book to <tt>quadl<\/tt>. And I think it never became very popular.\r\n<h4>Gauss, Kronrod, and <tt>quadgk<\/tt><a name=\"07a23423-d0df-4c19-a51e-3c163ad0b85f\"><\/a><\/h4>\r\nIn the mathematical software business when we were developing libraries using Fortran or C, we evaluated the cost of numerical quadrature by counting the number of function evaluations. In 2006 Larry Shampine, the SMU Professor who is a consultant to MathWorks and principal author of our ODE Suite, published the paper referenced below where he changed the ground rules for quadrature functions in the MATLAB environment. Since MATLAB function evaluations are vectorized, we should count the number of vectorized evaluations, not the number of individual evaluations. This is a better predictor of execution time.\r\n\r\nShampine's paper includes a MATLAB function that he calls <tt>quadva<\/tt>, where <tt>va<\/tt> stands for \"vectorized adaptive\". His code uses a 15-point Gauss-Kronrod formula. The Gauss formula is not only higher order than the Lobatto formula with the same number of points, it has the additional advantage that it does not require evaluating the integrand at the end points of an interval. This makes it possible to treat infinite intervals and singularites at the end points.\r\n\r\nWhen we added Shampine's <tt>quadva<\/tt> to MATLAB, we called it <tt>quadgk<\/tt>, for \"Gauss-Kronrod\".\r\n<h4>Comparison on <tt>humps<\/tt>.<a name=\"8f394fe9-1306-4888-ba9c-a2d792dad285\"><\/a><\/h4>\r\nLet's instrument <tt>humps<\/tt> to investigate the behavior of our quadrature routines on just one simple test case. The functions <tt>quad<\/tt> and <tt>quadl<\/tt> call the integrand once with a vector argument to initialize the integration and then call it more times with a shorter vector argument during the adaptive process. This is such an easy problem that <tt>quadgk<\/tt> gets all the information that it needs during the initialization.\r\n\r\nIn the following table, \"init\" is the length of the vector in the initialization call, \"kount\" is the number of times that <tt>humps<\/tt> is called after the initialization, \"length\" is the length of the vector argument in those calls, \"err\" is the error in the final result, and \"time\" is the execution time in milliseconds.\r\n\r\n$$\\begin{array}{lrrrrr}\r\n&amp; \\mbox{init} &amp; \\mbox{kount} &amp; \\mbox{length} &amp; \\mbox{err\\ \\ \\ \\ } &amp; \\mbox{time}\\\\\r\n\\mbox{quad} &amp; 7 &amp; 69 &amp; 2 &amp; 7.3 \\cdot 10^{-7} &amp; 0.87 \\\\\r\n\\mbox{quadl} &amp; 13 &amp; 37 &amp; 5 &amp; 1.8 \\cdot 10^{-10} &amp; 0.78 \\\\\r\n\\mbox{quadgk} &amp; 150 &amp; 0 &amp; - &amp; 2.5 \\cdot 10^{-14} &amp; 0.34\r\n\\end{array}$$\r\n\r\nWe see that the three functions perform very differently in this one experiment. The oldest function, <tt>quad<\/tt>, actually requires the fewest total function evaluations, if we count them in the traditional scalar sense, but it is the least accurate and the slowest. The newest function, <tt>quadgk<\/tt>, requires only one vectorized function evaluation, and is the fastest and the most accurate.\r\n<h4><tt>quadgk<\/tt> becomes <tt>integral<\/tt><a name=\"b20b2fd4-8c24-48a3-bfb8-804976454356\"><\/a><\/h4>\r\nIn 2012, with release 2012a, we came to the realization that most new MATLAB users who wanted to evaluate an integral had a hard time finding the function named \"quadgk\" that does the job. For all these years, we numerical analysts have been explaining to the rest of the world our quaint term <i>quadrature<\/i>. If that weren't enough, we have to tack on initials of mathematicians. So, we introduced <tt>integral<\/tt>. Actually <tt>quadgk<\/tt> is still available. <tt>Integral<\/tt> is just an easier to find and easier to use version of <tt>quadgk<\/tt>.\r\n<h4>Reference<a name=\"dc1fbbc4-e550-4971-9509-1fb0078e1939\"><\/a><\/h4>\r\nL.F. Shampine, \"Vectorized Adaptive Quadrature in MATLAB, <i>Journal of Computational and Applied Mathematics<\/i> 211 (2006), pp.131\u2013140, available at <a href=\"http:\/\/faculty.smu.edu\/shampine\/rev2vadapt.pdf\">&lt;http:\/\/faculty.smu.edu\/shampine\/rev2vadapt.pdf<\/a>&gt;\r\n\r\n<script>\/\/ <![CDATA[\r\nfunction grabCode_fd6db297af50441a9c61679cda6d4b9d() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='fd6db297af50441a9c61679cda6d4b9d ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' fd6db297af50441a9c61679cda6d4b9d';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2016 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('\r\n\r\n<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\r\n\r\n\r\n\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }\r\n\/\/ ]]><\/script>\r\n<p style=\"text-align: right; font-size: xx-small; font-weight: lighter; font-style: italic; color: gray;\">\r\n<a><span style=\"font-size: x-small; font-style: italic;\">Get\r\nthe MATLAB code<noscript>(requires JavaScript)<\/noscript><\/span><\/a>\r\n\r\nPublished with MATLAB\u00ae R2016a<\/p>\r\n\r\n<\/div>\r\n<!--\r\nfd6db297af50441a9c61679cda6d4b9d ##### SOURCE BEGIN #####\r\n%% Modernization of Numerical Integration, From Quad to Integral\r\n% The MATLAB functions for the numerical evaluation of integrals has\r\n% evolved from |quad|, through |quadl| and |quadgk|, to today's |integral|.\r\n\r\n%% Quadrature\r\n% Several years ago I regarded myself as an expert on numerical methods\r\n% for approximating integrals.  But then I stopped paying attention for a\r\n% while. So this is an opportunity to get reacquainted.  The task involves\r\n% a real-valued function $f(x)$ of a real variable $x$, defined on an\r\n% interval $a \\le x \\le b$. We seek to compute the value of the definite\r\n% integral,\r\n%\r\n% $$ \\int_a^b \\ f(x) \\ dx $$\r\n%\r\n\r\n%%\r\n% I wrote a chapter on the subject in my textbook _Numerical Computing\r\n% with MATLAB_.  Here is <https:\/\/www.mathworks.com\/moler.htmlquad.pdf a link % to that chapter>.  The chapter is titled \"Quadrature\".  One of the\r\n% first things I have to do in that chapter is explain that _quadrature_\r\n% is an old-fashioned term derived from the notion that the value of a\r\n% definite integral can be approximated by plotting the function on\r\n% graph paper and then counting the number of little squares that fall\r\n% underneath the graph.  Here is a coarse picture.\r\n%\r\n% <<quads.png>>\r\n%\r\n\r\n%% Adaptive Simpson's method and |quad|\r\n% For most of the time since the beginning of MATLAB, up until when I\r\n% wrote the chapter on quadrature twelve years ago, the MATLAB function for\r\n% evaluating definite integrals was one named |quad|.  It was based on\r\n% adaptive Simpson's method, pioneered by my grad school buddy\r\n% Bill McKeeman to show off recursion in Algol.\r\n\r\n%%\r\n% These two figures illustrate the basis for adaptive Simpson's method.\r\n%\r\n% <<simpsons.png>>\r\n%\r\n% The algorithm works recursively on subintervals of the original\r\n% interval [a,b].  If the values obtained from the 3-point Simpson's rule\r\n% and the 5-point composite Simpson's rule agree to within a specified\r\n% tolerance (which they clearly do not for these figures), then an\r\n% extrapolation taken from the two values is accepted as the value for\r\n% the integral over the subinterval.  If not, the algorithm is recursively\r\n% applied to the two halves of the subinterval.\r\n\r\n%% |Humps| and |quadgui|\r\n% I have used a function named |humps| to demonstrate numerical\r\n% integration for a long time.  It has a strong peak at x = 0.3 and a\r\n% milder peak  at x = 0.9.  It has found its way into the MATLAB |demos|\r\n% directory, so you have it on your own machine. A graph appears below.\r\n%\r\n% $$ h(x) = \\frac{1}{(x-0.3)^2+.01} + \\frac{1}{(x-0.9)^2+.04} - 6 $$\r\n%\r\n\r\n%%\r\n% An interactive graphical app named |quadgui| is included with the\r\n% programs for _Numerical Computing with MATLAB_, available from\r\n% <https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/37976-numerical-computing-with-matlab % MATLAB Central File Exchange>.  It shows |quad| in action.  Here is a\r\n% snapshot of |quad| part way through computing the integral of |humps|\r\n% over the interval [0,1].  Integration over the interval to the left of\r\n% 1\/4 is finished.  The integration over the subinterval [1\/4,3\/8]\r\n% has just been successfully completed.  A few samples to the right\r\n% of 3\/8 were computed during the initialization.  The integrand has been\r\n% evaluated 19 times.\r\n%\r\n% <<step19.png>>\r\n\r\n%%\r\n% Here is the final result.  The default tolerance for |quad| is $10^{-6}$\r\n% and the value of the integral obtained to six significant digits is\r\n% 29.8583.  The integrand has been evaluated 93 times.  Most of those\r\n% evaluations took place near the peaks where the |humps| varies rapidly.\r\n% This indicates that the adaptivity works as expected.\r\n%\r\n% <<humps.png>>\r\n%\r\n\r\n%%\r\n% |Humps| is a rational function.  It is possible to use the Symbolic\r\n% Toolbox to obtain the indefinite and then the definite integral\r\n% analytically.\r\n\r\nsyms x\r\nh = 1\/((x-3\/10)^2+1\/100) + 1\/((x-9\/10)^2+4\/100) - 6\r\nI = int(h,x)\r\nD = subs(I,1)-subs(I,0)\r\nformat long\r\nQ = double(D)\r\n\r\n%%\r\n% This confirms the six figure value shown by |quadqui|\r\n\r\n%% Lobatto, Kronrod and |quadl|\r\n% About the time in 2004 that I was writing the _Quadrature_ chapter for\r\n% _NCM_, the |quadl| function showed up.  This uses the _Lobatto_ rule.\r\n%\r\n% $$ \\int_{-1}^1 \\ f(x) \\ dx \\ \\approx \\ w_1 f(-1) + w_2 f(-x_1)\r\n%    + w_2 f(x_1) + w_1 f(1) $$\r\n%\r\n% with $x_1 = 1\/\\sqrt{5}, \\ w_1 = 1\/6, \\ w_2 = 5\/6$.\r\n\r\n%%\r\n% |Quadl| also employs the associated _Kronrod_ rule, which provides\r\n% 13 additional irrational abscissa and weights.  The Lobatto-Kronrod\r\n% rule is much higher order and consequently more accurate than Simpson's\r\n% rule.  As a result |quadl| is usually faster and more accurate than\r\n% |quad|.  But I just devoted a couple of exercises in the book to |quadl|.\r\n% And I think it never became very popular.\r\n\r\n%% Gauss, Kronrod, and |quadgk|\r\n% In the mathematical software business when we were developing libraries\r\n% using Fortran or C, we evaluated the cost of numerical quadrature by\r\n% counting the number of function evaluations.  In 2006 Larry Shampine,\r\n% the SMU Professor who is a consultant to MathWorks and principal author\r\n% of our ODE Suite, published the paper referenced below where he changed\r\n% the ground rules for quadrature functions in the MATLAB environment.\r\n% Since MATLAB function evaluations are vectorized, we should count the\r\n% number of vectorized evaluations, not the number of individual\r\n% evaluations.  This is a better predictor of execution time.\r\n\r\n%%\r\n% Shampine's paper includes a MATLAB function that he calls |quadva|,\r\n% where |va| stands for \"vectorized adaptive\".  His code uses a 15-point\r\n% Gauss-Kronrod formula.  The Gauss formula is not only higher order\r\n% than the Lobatto formula with the same number of points, it has the\r\n% additional advantage that it does not require evaluating the integrand\r\n% at the end points of an interval.  This makes it possible to treat\r\n% infinite intervals and singularites at the end points.\r\n\r\n%%\r\n% When we added Shampine's |quadva| to MATLAB, we called it |quadgk|,\r\n% for \"Gauss-Kronrod\".\r\n\r\n%% Comparison on |humps|.\r\n% Let's instrument |humps| to investigate the behavior of our quadrature\r\n% routines on just one simple test case.  The functions |quad| and |quadl|\r\n% call the integrand once with a vector argument to initialize the\r\n% integration and then call it more times with a shorter vector argument\r\n% during the adaptive process.  This is such an easy problem that |quadgk|\r\n% gets all the information that it needs during the initialization.\r\n\r\n%%\r\n% In the following table, \"init\" is the length of the vector in the\r\n% initialization call, \"kount\" is the number of times that |humps| is\r\n% called after the initialization, \"length\" is the length of the vector\r\n% argument in those calls, \"err\" is the error in the final result, and\r\n% \"time\" is the execution time in milliseconds.\r\n%\r\n% $$\\begin{array}{lrrrrr}\r\n%    & \\mbox{init}  & \\mbox{kount} & \\mbox{length} & \\mbox{err\\ \\ \\ \\ } & \\mbox{time}\\\\\r\n%    \\mbox{quad}    &   7 & 69     & 2 & 7.3 \\cdot 10^{-7} & 0.87 \\\\\r\n%    \\mbox{quadl}   &  13 & 37     & 5 & 1.8 \\cdot 10^{-10} & 0.78 \\\\\r\n%    \\mbox{quadgk}  & 150 &  0     & - & 2.5 \\cdot 10^{-14} & 0.34\r\n%   \\end{array}$$\r\n\r\n%%\r\n% We see that the three functions perform very differently in this one\r\n% experiment.  The oldest function, |quad|, actually requires the fewest\r\n% total function evaluations, if we count them in the traditional scalar\r\n% sense, but it is the least accurate and the slowest.  The newest\r\n% function, |quadgk|, requires only one vectorized function evaluation,\r\n% and is the fastest and the most accurate.\r\n\r\n%% |quadgk| becomes |integral|\r\n% In 2012, with release 2012a, we came to the realization that most\r\n% new MATLAB users who wanted to evaluate an integral had a hard time\r\n% finding the function named \"quadgk\" that does the job.  For all these\r\n% years, we numerical analysts have been explaining to the rest of the\r\n% world our quaint term _quadrature_.  If that weren't enough, we have to\r\n% tack on initials of mathematicians.  So, we introduced |integral|.\r\n% Actually |quadgk| is still available.  |Integral| is just\r\n% an easier to find and easier to use version of |quadgk|.\r\n\r\n%% Reference\r\n% L.F. Shampine, \"Vectorized Adaptive Quadrature in MATLAB,\r\n% _Journal of Computational and Applied Mathematics_ 211 (2006),\r\n% pp.131\u00e2\u20ac\u201c140, available at <http:\/\/faculty.smu.edu\/shampine\/rev2vadapt.pdf % http:\/\/faculty.smu.edu\/shampine\/rev2vadapt.pdf>\r\n\r\n##### SOURCE END ##### fd6db297af50441a9c61679cda6d4b9d\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/quads.png\" onError=\"this.style.display ='none';\" \/><\/div><!--introduction-->The MATLAB functions for the numerical evaluation of integrals has evolved from <tt>quad<\/tt>, through <tt>quadl<\/tt> and <tt>quadgk<\/tt>, to today's <tt>integral<\/tt>.\r\n\r\n<!--\/introduction-->... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2016\/05\/23\/modernization-of-numerical-integration-from-quad-to-integral\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[11,12,4,16,20],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/1576"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=1576"}],"version-history":[{"count":2,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/1576\/revisions"}],"predecessor-version":[{"id":2195,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/1576\/revisions\/2195"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=1576"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=1576"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=1576"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}