<< . .

. 59
( : 95)



. . >>

(PM) (TR) (CS) (RO) (AD)
m 4000 5600 250 50
3.24e ’ 10 3.30e ’ 10 2.98e ’ 10 3.58e ’ 11 3.18e ’ 10
Err
20007 29013 2519 5772 3540
flops

1/β
1 ’ K 2 x2 dx, with
TABLE 9.12. Methods for approximating I(f ) = 4π±
√ 0
±2 = (3 ’ 2 2)/100, β = 10 and K 2 = β 2 (1 ’ ±2 β 2 )




5
10
1


0.8

4
10
0.6


0.4

3
10
0.2


0
2
10
0.1
0.06 0.08
0.04
0 0.02
’0.2
0.02 0.04 0.06 0.1
0.08
0


FIGURE 9.7. Distribution of quadrature nodes (left); integration stepsize density
in the approximation of integral (9.62) (right)




9.10.2 Computation of the Wind Action on a Sailboat Mast
Let us consider the sailboat schematically drawn in Figure 9.8 (left) and
subject to the action of the wind force. The mast, of length L, is denoted by
the straight line AB, while one of the two shrouds (strings for the side sti¬-
ening of the mast) is represented by the straight line BO. Any in¬nitesimal
element of the sail transmits to the corresponding element of length dx of
the mast a force of magnitude equal to f (x)dx. The change of f along with
the height x, measured from the point A (basis of the mast), is expressed
by the following law
±x ’γx
f (x) = e ,
x+β

where ±, β and γ are given constants.
9.10 Applications 411

The resultant R of the force f is de¬ned as
L

f (x)dx ≡ I(f ),
R= (9.63)
0

and is applied at a point at distance equal to b (to be determined) from
the basis of the mast.

B
B
mast
shroud


wind

1111111111
0000000000
direction

0000000000
1111111111 L
R
1111111111
0000000000 f dx
T
dx
T
b
M
A
O
O A
H
V

FIGURE 9.8. Schematic representation of a sailboat (left); forces acting on the
mast (right)

Computing R and the distance b, given by b = I(xf )/I(f ), is crucial for the
structural design of the mast and shroud sections. Indeed, once the values
of R and b are known, it is possible to analyze the hyperstatic structure
mast-shroud (using for instance the method of forces), thus allowing for the
computation of the reactions V , H and M at the basis of the mast and the
traction T that is transmitted by the shroud, and are drawn in Figure 9.8
(right). Then, the internal actions in the structure can be found, as well as
the maximum stresses arising in the mast AB and in the shroud BO, from
which, assuming that the safety veri¬cations are satis¬ed, one can ¬nally
design the geometrical parameters of the sections of AB and BO.
For the approximate computation of R we have considered the compos-
ite midpoint, trapezoidal and Cavalieri-Simpson rules, denoted henceforth
by (MP), (TR) and (CS), and, for a comparison, the adaptive Cavalieri-
Simpson quadrature formula introduced in Section 9.7.2 and denoted by
(AD). Since a closed-form expression for the integral (9.63) is not available,
the composite rules have been applied taking mk = 2k uniform partitions
of [0, L], with k = 0, . . . , 15.
We have assumed in the numerical experiments ± = 50, β = 5/3 and
γ = 1/4 and we have run Program 77 taking tol=10’4 and hmin=10’3 .
The sequence of integrals computed using the composite formulae has been
stopped at k = 12 (corresponding to mk = 212 = 4096) since the remaining
412 9. Numerical Integration

0
10
’1
10
’2
10
’3
(TR)
10
’4
10 (PM)
’5
10
’6
10
(CS)
’7
10
’8
10 (AD)
’9
10
0 20 40 60 80 100 120


FIGURE 9.9. Relative errors in the approximate computation of the integral
(±xe’γx )/(x + β)dx
L
0



elements, in the case of formula CS, di¬er among them only up to the last
signi¬cant ¬gure. Therefore, we have assumed as the exact value of I(f )
(CS)
the outcome I12 = 100.0613683179612 provided by formula CS.
(CS)
We report in Figure 9.9 the logarithmic plots of the relative error |I12 ’
Ik |/I12 , for k = 0, . . . , 7, Ik being the generic element of the sequence for
the three considered formulae. As a comparison, we also display the graph
of the relative error in the case of formula AD, applied on a number of
(nonuniform) partitions equivalent to that of the composite rules.
Notice how, for the same number of partitions, formula AD is more accu-
rate, with a relative error of 2.06 · 10’7 obtained using 37 (nonuniform)
partitions of [0, L]. Methods PM and TR achieve a comparable accuracy
employing 2048 and 4096 uniform subintervals, respectively, while formula
CS requires about 64 partitions. The e¬ectiveness of the adaptivity pro-
cedure is demonstrated by the plots in Figure 9.10, which show, together
with the graph of f , the distribution of the quadrature nodes (left) and the
function ∆h (x) (right) that expresses the (piecewise constant) density of
the integration stepsize h, de¬ned as the inverse of the value of h over each
active interval A (see Section 9.7.2).
Notice also the high value of ∆h at x = 0, where the derivatives of f are
maximum.



9.11 Exercises
1. Let E0 (f ) and E1 (f ) be the quadrature errors in (9.6) and (9.12). Prove
that |E1 (f )| 2|E0 (f )|.
2. Check that the error estimates for the midpoint, trapezoidal and Cavalieri-
Simpson formulae, given respectively by (9.6), (9.12) and (9.16), are special
instances of (9.19) or (9.20). In particular, show that M0 = 2/3, K1 =
9.11 Exercises 413

20 30


25
15


20
10

15

5
10


0
5


’5 0
0 2 4 6 8 10 0 2 4 6 8 10


FIGURE 9.10. Distribution of quadrature nodes (left); integration step density
in the approximation of the integral 0 (±xe’γx )/(x + β)dx (right)
L




’1/6 and M2 = ’4/15 and determine, using the de¬nition, the degree of
exactness r of each formula.
b
[Hint: ¬nd r such that In (xk ) = a xk dx, for k = 0, . . . , r, and In (xj ) =
bj
x dx, for j > r.]
a
n
3. Let In (f ) = k=0 ±k f (xk ) be a Lagrange quadrature formula on n + 1
nodes. Compute the degree of exactness r of the formulae:
(a) I2 (f ) = (2/3)[2f (’1/2) ’ f (0) + 2f (1/2)],
(b) I4 (f ) = (1/4)[f (’1) + 3f (’1/3) + 3f (1/3) + f (1)].
Which is the order of in¬nitesimal p for (a) and (b)?
[Solution: r = 3 and p = 5 for both I2 (f ) and I4 (f ).]
4. Compute df [x0 , . . . , xn , x]/dx by checking (9.22).
[Hint: proceed by computing directly the derivative at x as an incremental
ratio, in the case where only one node x0 exists, then upgrade progressively
the order of the divided di¬erence.]

1
5. Let Iw (f ) = 0 w(x)f (x)dx with w(x) = x, and consider the quadrature
formula Q(f ) = af (x1 ). Find a and x1 in such a way that Q has maximum
degree of exactness r.
[Solution: a = 2/3, x1 = 3/5 and r = 1.]
6. Let us consider the quadrature formula Q(f ) = ±1 f (0) + ±2 f (1) + ±3 f (0)
1
for the approximation of I(f ) = 0 f (x)dx, where f ∈ C 1 ([0, 1]). Determine
the coe¬cients ±j , for j = 1, 2, 3 in such a way that Q has degree of
exactness r = 2.
[Solution: ±1 = 2/3, ±2 = 1/3 and ±3 = 1/6.]
7. Apply the midpoint, trapezoidal and Cavalieri-Simpson composite rules to
approximate the integral
1
|x|ex dx,
’1

and discuss their convergence as a function of the size H of the subintervals.
414 9. Numerical Integration

1
8. Consider the integral I(f ) = 0 ex dx and estimate the minimum number m
of subintervals that is needed for computing I(f ) up to an absolute error
¤ 5 · 10’4 using the composite trapezoidal (TR) and Cavalieri-Simpson
(CS) rules. Evaluate in both cases the absolute error Err that is actually
made.
[Solution: for TR, we have m = 17 and Err = 4.95 · 10’4 , while for CS,
m = 2 and Err = 3.70 · 10’5 .]
9. Consider the corrected trapezoidal formula (9.30) and check that |E1 (f )|
corr
corr
4|E2 (f )|, where E1 (f ) and E2 (f ) are de¬ned in (9.31) and (9.16), re-
spectively.
10. Compute, with an error less than 10’4 , the following integrals:

sin(x)/(1 + x4 )dx;
(a) 0

e’x (1 + x)’5 dx;
(b) 0
2

cos(x)e’x dx.
(c) ’∞

11. Use the reduction midpoint and trapezoidal formulae for computing the
y
dxdy over the domain „¦ = (0, 1)2 . Run
double integral I(f ) = „¦
(1 + xy)
Programs 78 and 79 with M = 2i , for i = 0, . . . , 10 and plot in log-scale
the absolute error in the two cases as a function of M . Which method is
the most accurate? How many functional evaluations are needed to get an
(absolute) accuracy of the order of 10’6 ?
[Solution: the exact integral is I(f ) = log(4) ’ 1, and almost 2002 = 40000
functional evaluations are needed.]
10
Orthogonal Polynomials in
Approximation Theory




Trigonometric polynomials, as well as other orthogonal polynomials like
Legendre™s and Chebyshev™s, are widely employed in approximation theory.
This chapter addresses the most relevant properties of orthogonal poly-
nomials, and introduces the transforms associated with them, in particular
the discrete Fourier transform and the FFT, but also the Zeta and Wavelet
transforms.
Application to interpolation, least-squares approximation, numerical dif-
ferentiation and Gaussian integration are addressed.




10.1 Approximation of Functions by Generalized
Fourier Series
Let w = w(x) be a weight function on the interval (’1, 1), i.e., a nonneg-
ative integrable function in (’1, 1). Let us denote by {pk , k = 0, 1, . . . } a
system of algebraic polynomials, with pk of degree equal to k for each k,
mutually orthogonal on the interval (’1, 1) with respect to w. This means
that

1

pk (x)pm (x)w(x)dx = 0 if k = m.
’1
416 10. Orthogonal Polynomials in Approximation Theory

1 1/2
Set (f, g)w = ’1 f (x)g(x)w(x)dx and f w = (f, f )w ; (·, ·)w and · w
are respectively the scalar product and the norm for the function space
1
f : (’1, 1) ’ R, f 2 (x)w(x)dx < ∞ .
L2 L2 (’1, 1)
= = (10.1)
w w
’1

For any function f ∈ L2 the series
w

+∞
(f, pk )w
Sf = fk p k , with fk = ,
pk 2 w
k=0


is called the generalized Fourier series of f , and fk is the k-th Fourier
coe¬cient. As is well-known, Sf converges in average (or in the sense of
L2 ) to f . This means that, letting for any integer n
w

n
fn (x) = fk pk (x) (10.2)
k=0

(fn ∈ Pn is the truncation of order n of the generalized Fourier series of
f ), the following convergence result holds

f ’ fn
lim = 0.
w
n’+∞

Thanks to Parseval™s equality, we have
+∞
2 2 2
f = fk pk
w w
k=0

+∞
and, for any n, f ’fn 2 = k=n+1 fk pk 2 is the square of the remainder
2
w w
of the generalized Fourier series.
The polynomial fn ∈ Pn satis¬es the following minimization property

f ’ fn = min f ’ q w. (10.3)
w
q∈Pn

+∞
Indeed, since f ’ fn = k=n+1 fk pk , the property of orthogonality of
polynomials {pk } implies (f ’ fn , q)w = 0 ∀q ∈ Pn . Then, the Cauchy-
Schwarz inequality (8.29) yields

f ’ fn = (f ’ fn , f ’ fn )w = (f ’ fn , f ’ q)w + (f ’ fn , q ’ fn )w
2
w

= (f ’ fn , f ’ q)w ¤ f ’ fn f ’q ∀q ∈ Pn ,
w,
w

and (10.3) follows since q is arbitrary in Pn . In such a case, we say that fn
is the orthogonal projection of f over Pn in the sense of L2 . It is therefore
w
interesting to compute the coe¬cients fk of fn . As will be seen in later

<< . .

. 59
( : 95)



. . >>