Integrating a Function
Instead of finding the area under a curve defined by a set of data points, you may wish to integrate a function F(x). You could simply create a table of function values and use one of the methods described in earlier sections to calculate the area. But a more convenient solution would be to create a custom function that uses the Formula property of the cell to get the worksheet formula to be integrated, in the same way that was used in the preceding chapter, and uses the formula to find the area under the curve. This approach will be described in subsequent sections.
Integrating a Function
Defined by a Worksheet Formula by Means of a VBA Function Procedure
In this section, the trapezoidal and Simpson's rule methods are implemented as VBA custom functions, using an approach similar to that used in the differentiation functions of the previous chapter. The Formula property of the cell is used to get the worksheet formula to be differentiated into the VBA code as text. Then the SUBSTITUTE worksheet function is used to replace the variable of interest by an incremented value, and the Evaluate method used to get the new value of the formula. These values are used to calculate the area of each panel, and the areas of the panels are summed to obtain the area under the curve.
This function procedure can be used to integrate an expression F(x) defined by a worksheet formula, between specified lower and upper limits a and b respectively. A table of function values is not required.
The syntax of the function is In teg rate (express/on, variable, fromjower, to_upper). The argument expression is the integrand, the expression to be integrated. The argument variable is the variable of integration, and the arguments fromjower and to_upper are the lower and upper limits of integration, respectively. The VBA code is shown in Figure 78. Function procedures for both trapezoidal (IntegrateT) and Simpson's rule (Integrates) methods are shown.
The range of x values over which the integration is to be performed (to_upper  fromjower) is divided into N panels. The user can adjust the accuracy of the integration by changing the value of N in the procedure, with a concomitant increase in calculation time.
Option Explicit
Function lntegrateT(expression, variable, fromjower, to_upper) 'Simple trapezoidal area integration
Dim FormulaString As String, T As String, Xref As String
Dim H As Double, area As Double, X As Double
Dim N As Integer, K As Integer, J As Integer
Dim NRepI As Integer
Dim temp As String
Dim F1 As Double, F2 As Double
FormulaString = expression.Formula
T = Application.ConvertFormula(FormulaString, xlA1, xlA1, xIAbsolute) XRef = variable. Address
H = (to_upper  fromjower) / N area = 0 X = fromjower
NRepI = (Len(T)  Len(Application.Substitute(T, XRef,""))) / Len(XRef) For K = 1 To N
For J = NRepI To 1 Step 1 temp = Application.Substitute(T, XRef, X & "", J) If lsError(Evaluate(temp)) Then GoTo pt1 T = temp pt1: Next J F1 = Evaluate (T)
T = Application.ConvertFormula(FormulaString, xlA1, xlA1, xlAbsolute) For J = NRepI To 1 Step 1 temp = Application.Substitute(T, XRef, X + H & "", J) If lsError(Evaluate(temp)) Then GoTo pt2 T = temp pt2: Next J F2 = Evaluate(T)
IntegrateT = area
End Function_
Figure 78. VBA Function procedure to integrate a worksheet formula by the trapezoidal approximation method, (folder 'Chapter 07 Examples,' workbook 'Integration,' module 'Simplelntegration')
Function lntegrateS(expression, variable, fromjower, to_upper) 'Simpson's 1/3 rule area integration
Dim FormulaString As String, T As String, Xref As String Dim H As Double, area As Double, X As Double Dim N As Integer, K As Integer, J As Integer Dim NRepI As Integer Dim temp As String
Dim Y0 As Double, Y1 As Double, Y2 As Double
FormulaString = expression.Formula XRef = variable.Address
T = Application.ConvertFormula(FormulaString, xlA1, xlA1, xIAbsolute) NRepI = (Len(T)  Len(Application.Substitute(T, XRef,""))) / Len(XRef) For J = NRepI To 1 Step 1 temp = Application.Substitute(T, XRef, fromjower + X & "", J) If lsError(Evaluate(temp)) Then GoTo pt1 T = temp pt1: Next J
T = Application.ConvertFormula(FormulaString, xlA1, xlA1, xIAbsolute) For J = NRepI To 1 Step 1 temp = Application.Substitute^, XRef, fromjower + X + H & "", J) If lsError(Evaluate(temp)) Then GoTo pt2 T = temp pt2: Next J Y1 = Evaluate (T)
T = Application.ConvertFormula(FormulaString, xlA1, xlA1, xIAbsolute) For J = NRepI To 1 Step 1 temp = Application.Substitute(T, XRef, fromjower + X + 2 * H & "", J) If lsError(Evaluate(temp)) Then GoTo pt3 T = temp pt3: Next J Y2 = Evaluate (T)
Integrates = area
End Function_
Figure 79. VBA function procedure to integrate a worksheet formula by Simpson's method, (folder 'Chapter 07 Examples', workbook 'Integration', module 'Simplelntegration')
Some results returned by the IntegrateT and Integrates functions are shown in Figures 710 and 711, respectively. In general, results are more accurate when using the Simpson's method function.
A 
B 
C 
D 
E 
F G 
^ H I';]  
1 
Integration 
function using simple trapezoidal appro*.  
2 
Function 
X 
Fix) 
from 
to 
Area exact 
% error 
3 
X3 
1 
1 
0 
1 
0.2500003 0.25 
1.0E04 
4 
4<0r) 
1 
0 
0 
1 
3.14156 3.14159 
1 2E03 
5 
triangle' 
1 
2 
0 
2 
4.00000000 4 
O.OE+OO 
6 
Gaussian 
95 
0.035 
20 
180 
1.00000000 1 
3.1E12 
7  
8 
* 
slope 
= 2, intercept= 4  
9 
** ¡I 
= 100 
,o= 10 
A 
8 
C 
D 
E 
F  G 
H  
1 
Integration function using Simpson's method  
2 
Function 
X 
FOO 
from 
to 
Area 
exact 
% error 
3 
x3 
1 
1 
0 
1 
0,2500000 
0,25 
4.4E1 4 
4 
4iCIx2) 
1 
0 
0 
1 
3.141588 
3.141593 
1.6E04 
5 
triangle 
1 
2 
0 
2 
4.0000000 
4 
2.2E14 
6 
Gaussian 
95 
0.035 
20 
180 
1.0000000 
1 
1.3E13 
* '  
8 
* 
slope 
= 2, intercept = 4  
9 
** ^ 
= 100 
, u= 10 
Figure 711. Some results returned by the Integrates custom function, (folder 'Chapter 07 Examples', workbook 'Integration', sheet Simpson Integration Fn')
Figure 711. Some results returned by the Integrates custom function, (folder 'Chapter 07 Examples', workbook 'Integration', sheet Simpson Integration Fn')
Because some functions may require a large number of iterations, there may be a noticeable delay in calculation.
Gaussian Quadrature
The preceding methods for numerical integration employ evenly spaced values of x at which the function is evaluated. Other formulas have been developed whereby the function is evaluated at specially selected values of x. These Gaussian quadrature formulas are significantly more efficient, in terms of the accuracy of the evaluation.
Gaussian quadrature formulas involve the evaluation of the function at a set of Xj values (nodes), with the use of a corresponding set of weights w„ in the following formula
The nodes and weights can be derived from certain kinds of polynomials. The Legendre polynomials will be used here to determine the values of x{ and Wj. The Legendre polynomials are a set of polynomials of degree N. Increasing N provides an increase in accuracy of evaluation but requires a concomitant increase in computation time. Values of Legendre polynomials for N up to 100 have been published.
The integration need not be limited solely to the interval 1 to 1. By employing a change of variable
the integral expression is b , 1 , A= JF(x)dx =j>
and equation 79 becomes b ba^ J (b  a)zi +(b + a)
a which permits integration over any range.
The code shown in Figure 712 performs Gaussian quadrature using equation 712 and a tenthorder Legendre polynomial. Some results returned by the function are shown in Figure 713.
Option Explicit
Function lntegrate(expression, variable, fromjower, to_upper, Optional _ tolerance)
Dim FormulaString As String, XAddress As String Dim result As Double
FormulaString = expression.Formula XAddress = variable.Address 'Default is absolute FormulaString = Application.ConvertFormula(FormulaString, xlA1, xlA1, _ xIAbsolute)
Call GaussLeGendre10(FormulaString, XAddress, fromjower, to_upper, _
tolerance, result) Integrate = result End Function
Sub GaussLeGendre10(expression, XRef, fromjower, to_upper, tolerance, result)
'Uses tenpoint GaussLegendre quadrature formula. 'Adapted from Shoup, p.203
Dim XJ As Variant, AJ As Variant
Dim TotalArea As Double, OldArea As Double, area As Double Dim T As String, temp As String
Dim I As Integer, J As Integer, K As Integer, JJ As Integer Dim N As Integer, NRepI As Integer
Dim A As Double, B As Double, C As Double, D As Double, F As Double Dim H As Double
XJ = Array(0.973906528517172, 0.865063366688984, 0.679409568299024, 0.433395394129247, 0.148874338981631, 0.973906528517172,
0.865063366688984, 0.679409568299024, 0.433395394129247, 0.148874338981631) AJ = Array(0.066671344308688, 0.149451349915058, 0.219086362515982, 0.269266719309996, 0.295524224714753, 0.066671344308688, 0.149451349915058, 0.219086362515982, 0.269266719309996, 0.295524224714753)
If IsMissing(tolerance) Then tolerance = 0.0000000001 OldArea = 0 N = 1
For K = 1 To 10 'increments divided by 1,2,4,8,16,32,64,128,256,512 area = 0
A = from Jower + (I 1) • H B = A + H C = (B + A) / 2 D = (B  A) / 2
NRepI = (Len(T)  Len(Application.Substitute(T, XRef,""))) / Len(XRef) For J J = NRepI To 1 Step1 temp = Application.Substitute^", XRef, C + D * XJ(J) & "", JJ) If lsError( Evaluate(temp)) Then GoTo pt1 T = temp pt1: Next JJ F = Evaluate(T) area = area + AJ(J) * F Next J Next I
If Abs((area  OldArea) / area) < tolerance Then GoTo AIIDone
OldArea = area
Next K
AIIDone:
result = area
End Sub
Figure 712. Integrate custom function, (folder 'Chapter 07 Examples', workbook 'Integration', module 'Legendreintegration')
A 
B I 
c 
D 
E 
F i G 
.H  
1 
Custom Function to Integrate Aiea Under a Curve  
2 
Integ rate {formu5a_in_ceS5, variable, from_ lower; fo_ upper, tolerance)  
3 
Function 
X 
FOO 
ftom 
to 
area exact value 
% error 
4 
4V(1 >r)dx 
1 
0 
1 
1 
3.1415928 3.1415927 
4.0E06 
S 
x\jx 
1 
1 
0 
1 
0.2500000002 0,25 
7.SE08 
6 
(mx+b)dx 
2 
0 
0 
2 
4.000000003 4 
7.6E08 
7 
sin(x)dx 
o 
0 
0 
3.14 
2.000000001 2 
7.4E08 
9 
(x e"')dx 
1 
0.37 
0 
1000 
1000000001 1 
9.5E08 
9 
(1A0 dx 
1 
1 
1 
2 
0.69314718 0.69314718 
7.7E08 
Figure 713. Some results returned by the Integrate custom function, (folder 'Chapter 07 Examples', workbook 'Integration', sheet 'GaussLegendre Integration Fn')
Figure 713. Some results returned by the Integrate custom function, (folder 'Chapter 07 Examples', workbook 'Integration', sheet 'GaussLegendre Integration Fn')
Early versions of this program returned inaccurate results when the range h — a was large. The function Integrate illustrates one approach to overcoming this problem. First, the integral is evaluated over the total range b  a. Then the interval is divided into two halves and each "panel" is integrated separately. The sum of the two panels is compared to the previous value. If the difference is larger than a tolerance value, the interval is divided into quarters, the areas summed and so on. The process is continued for 10 cycles of iteration (512 panels) or until the area difference is less than a specified tolerance.
Because some functions may require a large number of iterations, there may be a noticeable delay in calculation. Increasing the value of tolerance should speed up calculation, but only at the expense of accuracy.
Integration with an Upper or Lower Limit of Infinity
Integrals such as
a can be evaluated by summing the areas of a number of panels covering the range from x = a to x = a suitably large value. It is to be expected that as x—> oo the area of panel(x) » zero. Thus the integral can be evaluated by summing the integrals of a series of panels of increasing width (e.g., from 01, 110, 10100, etc), ending the summation when the area of the last panel is suitably small. Manual adjustment of the panel widths is easily done by inspection of the results. Figure 714 shows a typical result.
1 A „ 
B 
C 
D 
E 
F i 
G 
H !  
32  
Integrating to upper limit of infinity  
33 
Function 
X 
m 
from 
to 
integrated value 
I 
% error  
34 
y = 1 ¿£(1 +x)*sqtt(x) 
1 0.5 
0 
0.01 
0.1989708 : 
0.1989708 
93.67  
35 
1 
0.5 
0.01 
01 
0.4132174 
0.6121882 
80.51  
36 
1 
0.5 
IM: 
1 
0.9582416 
1.5704298 
50.01  
37 
1 
0.5 
1 
10 
0.9582416 
2.5286714 
19 51  
36 
1 I 0.5 
10 
100 
0.4132174 
2.9418888 
6.36  
39 
1 
0.5 
100 
1000 
0.1361128 
3.0780016 
2.02  
40 
0.5 
1E+03 
1E+04 
0.0432252 
3.1212268 
0.65  
41 
1 
0.5 
1E+04 
1E+05 
0.0136748 
3.1349016 
0.21  
42 
1 
0.5 
1E+05 
1E+06 
0.0043245 
3.1392261 
0.08  
43~ 
1 
0.5 
1E+06 
1E+08 
0.0018000 
3.1410261 
0.02  
44 
1 
0.5 
1E+08 
1 E+10 
0.0001800 
3.1412061 
0.01  
45 
1 
0.5 
1E+10 
1E+14 
0,0000198 
3.1412259 
0.01  
46 
1 
0.5 
1E+14 
1E+18 
0,0000002 
3.1412261 
0.01  
47 
exact value 
3.1415927 
Figure 714. Integrating from a lower limit to an upper limit of infinity. Results returned by the Integrate custom function, (folder 'Chapter 07 Examples', workbook 'Integration', sheet 'Integrating to infinity by sum')
Distance Traveled Along a Curved Path
The length of a plane curve can be estimated by dividing the curve into segments, as in Figure 715, and approximating the length of the curve segment by the straight line AB. The length of AB = V(Ax)2 + (Ay)2 . The distance along
A 
B _ 
c  
1 
Distance Travelled Along a Curve  
2 
1/4 circle with r = 
= 1  
3 
(distance should be n/2)  
4 
X 
y 
d 
5 ' 
0.000 i 1.000  
6 
0.050 
0.999 
0.050 
7 ' 
0 075 
0.997 
0.G25 
8 
0 100 
0.995 
0.025 
9 
0.125 
0.992 
0.025 
51 
0.980 
0.199 
0.024 
52 
0.985 
0.173 
0.027 
53 
0.990 
0.141 
0.032 
54 
0.9925 
0.122 
0.019 
55 
0.9950 
0.100 
0.023 
56 
0.9975 
0.071 
0.029 
57 
0.999 
0.045 
0.026 
58 
1.000 
0.000 
0 045 
59 
I Sumx2 = 
3.14145  
60 
% error = 
4 JE03 
Figure 716. Approximating the circumference of a circle of radius 1.
Note that the rows between 9 and 51 are hidden, (folder 'Chapter 07 Examples', workbook 'Curve Distance', sheet 'Curve Distance (Circle)')
Figure 716. Approximating the circumference of a circle of radius 1.
Note that the rows between 9 and 51 are hidden, (folder 'Chapter 07 Examples', workbook 'Curve Distance', sheet 'Curve Distance (Circle)')
The procedure is illustrated by estimating the length of one quarter of a circle of radius r= 1. The equation of the circle is x2 + y2 = 1, or y =%\x2 . As shown in Figure 716, the value of y and the distance d between successive points was calculated from x = 0 to * = 1, using an x increment of 0.025. Near the end of the range of x values, where y changes more rapidly, the x increment was decreased. The formula in cell C6 is
=SQRT((A8A5)A2+(B8B5)A2) The sum of the distances x 2, in cell C59 is a reasonable estimate of %.
Problems
Answers to the following problems are found in the folder "Ch. 07 (Integration)" in the "Problems & Solutions" folder on the CD.
00 x2
method.
2. Integrate the following expressions, using one of the custom functions for integration.
o lnx
3. Evaluate the elliptic integral nil j\/l  (1/2)sin2 xdx
4. An ellipse is a plane figure described by the locus of a point P(x, y) that moves such that the sum of its distances from two fixed points (foci) is a constant. If the ellipse has foci located at A (c, 0) and B (c, 0) and the distance ACB is 2a, then by setting b = Va2 c2 , the equation of the ellipse is simplified to
(a and b are termed the semiaxes of the ellipse).
Figure 717. Approximating the circumference of an ellipse.
For the ellipse shown in Figure 717, with foci at x = 0.5, y = 0 and x = 0.5, y = 0 and a= 1, determine the circumference of the ellipse.
5. Determine the area of the ellipse of problem 74.
6. Find the area between the curve y  2x  x2 and the line y = 3.
7. Find the area between the curve y = 2x  xz and the line y = 2.5x 2.3.
8. Find the area enclosed between the two curves shown in Figure 77: y\ = x3 
20x2  IOOjc + 2000 and_y2 = 2x3  5x2  300x + 1000. The curves intersect in the region between x = 5 and x = 15.
(a and b are termed the semiaxes of the ellipse).
Figure 717. Approximating the circumference of an ellipse.
9. The area between the curve y = x2 and the horizontal line y = 4 is divided into two equal areas by the horizontal line y — c. Find c.
10. The area between the curve y = x2 + 3 and the line y = 12 is divided into two equal areas by the line y — c. Findc.
11. Integrate the following expression.
12. Integrate the following expressions, using the custom function for integration.
00 j
This Page Intentionally Left Blank
Responses

Vittorio6 years ago
 Reply

phillipp bader5 years ago
 Reply

poppy5 years ago
 Reply

vihtori5 years ago
 Reply

lalli4 years ago
 Reply

zewdi asmara4 years ago
 Reply