Il problema della quadratura
 Data una funzione f(x) definita in un intervallo [a,b], si
vuole valutare l’integrale:
b
I   f ( x)dx
a
a partire dai valori della funzione integranda f(x) in un
insieme di punti compresi nell’intervallo di integrazione
 In tutti i metodi di quadratura si effettua una somma di
valori della funzione integranda
 Un buon metodo di quadratura deve:
 valutare l’integrale con la maggior precisione possibile
 sfruttare il minor numero possibile di valori della funzione
integranda
Notazioni
 Supponiamo di avere una serie di N ascisse equispaziate x1, x2, ... , xN :
 x1=a; xN=b
 h = distanza tra ciascuna coppia di ascisse
 (b-a) = (N-1)h
 xk=xk+(k-1)h con k=1,...,N
 poniamo f(xi)=fi
 Formule chiuse:
 utilizzano nel calcolo i valori di f1 e fN
 Formule aperte:
 non utilizzano nel calcolo uno o entrambi i valori di f1 e fN
 possono essere utili se il valore di f in uno degli estremi di integrazione è infinito
(purché la singolarità sia integrabile)
y
fi=f(xi)
h
x1=a
x2
xi
xN=b
x
Regola del trapezio
 Consiste nell’approssimare l’integrale nell’intervallo tra xj e xj+1
nel modo seguente:
x j 1

xj

1
1

f ( x)dx  h  f j  f j 1   O h3 f 
2
2


 Se f(x)≥0, tale valore rappresenta l’area del trapezio di basi f j e f j+1
e altezza h
 in sostanza, nell’intervallo tra xj e xj+1 la f(x) viene approssimata
da un polinomio di primo grado
 Il valore dell’integrale calcolato con la regola del trapezio
differisce dal valore vero per un termine che è dell’ordine di h3
per la derivata seconda della funzione calcolata in un punto
(non noto) dell’intervallo [xj,xj+1]
 La formula del trapezio è esatta per polinomi fino al primo
grado
Regola del trapezio estesa
 Utilizziamo la regola del trapezio N-1 volte negli intervalli
[x1,x2], [x2,x3], ..., [xN-1,xN]:
xN

x1


1  1
1 
1 
1
1
f ( x)dx  h  f1  f 2   h  f 2  f 3     h  f N 1  f N   N  1O h 3 f  
2  2
2 
2 
2
2
 b  a 3 
1 
1
h  f1  f 2  f 3    f N 1  f N    N  1O
f   
3
2 
2
  N  1

 b  a 3 
1 
1
h  f1  f 2  f 3    f N 1  f N   O
f  
2
2 
2
 N

 nel calcolo precedente si è assunto N-1≈N, il che è vero per N
grande
 La precisione migliora con il quadrato del numero di punti
utilizzati per il calcolo
 raddoppiando i punti l’errore diminuisce di un fattore 4
Applicazione della regola del trapezio
1
iterazioni
2
3
4
...
 Si procede per approssimazioni successive:
 nella prima iterazione si utilizzano i valori di f(x) negli estremi di
integrazione
 nella n-esima iterazione vengono aggiunti 2n-2 valori di f(x) in
corrispondenza dei punti medi degli intervalli elementari definiti
dall’iterazione precedente
 alla n-esima iterazione l’intervallo di integrazione risulta diviso in 2n-1
intervalli elementari
Implementazione dell’algoritmo (1)
 Poniamo b-a=Δ e indichiamo con In l’integrale calcolato
nella n-esima iterazione.
 Alla prima iterazione si ha:
1
1

I1    f a   f b 
2
2

 Alla seconda iterazione avremo:
I2 
 1
 1
 1
 






f
a

f
a


f
b

I

f
a





 2 1 2
2  2
2
2
2





 Alla terza iterazione si avrà:
I3 

 1


3  1







f
a

f
a


f
a


f
a


f
b







4  2
4
2
4
2







1
 

3 

I2   f  a    f  a 

2
4 
4
4


Implementazione dell’algoritmo (2)
 Generalizzando il risultato trovato in precedenza
possiamo concludere che:
n 2
1
 2
I n  I n 1  n 1 
2
2 i 1

 

f  a  n 1  i  1 n  2 
2
2 

 Poiché nel passare da un’iterazione alla successiva il
numero di punti utilizzati per il calcolo raddoppia, la
precisione del calcolo migliora di un fattore 4.
 In generale, la procedura iterativa viene fermata
quando:
I n  I n1   I n1
dove ε è la precisione richiesta per il calcolo
Esempi (1)
 Supponiamo di voler calcolare l’integrale
1

dx

 0.785398
2
1 x
4
0
1
I 
con una precisione ε=10-6. Si ha:
Iterazione (n)
In
1
0.75
2
0.775
3
0.782794
4
0.784747
5
0.785235
6
0.785357
7
0.785388
8
0.785396
9
0.785398
10
0.785398
Esempi (2)
 Supponiamo di voler calcolare l’integrale
2
1
I   dx  ln 2  0.693147 
x
1
con una precisione ε=10-6. Si ha:
Iterazione (n)
In
1
0.75
2
0.708333
3
0.697024
4
0.694122
5
0.693391
6
0.693208
7
0.693162
8
0.693151
9
0.693148
10
0.693147
11
0.693147
Regola di Simpson (1)
 Consideriamo l’intervallo tra xj-1 e xj+1
 Usiamo la seguente notazione:
 xj-1=xj-h; f(xj-1)=f j-1
 f(xj)=f j
 xj+1=xj+h; f(xj+1)=f j+1
 Nell’intervallo in esame approssimiamo la f(x) con una
parabola passante per i tre punti (xj-1,fj-1), (xj,fj), (xj+1,fj+1):
f ( x)  f j
 ax  x   bx  x 
2
j
j
 scrivendo l’equazione della parabola in questa forma è
automaticamente rispettata la condizione f(xj)=f j
 I coefficienti a e b si determinano imponendo le condizioni;
 f(xj-h)=f j-1
 f(xj+h)=f j+1
Regola di Simpson (2)
y
fj
fj+1
f j-1
h
xj-1
h
xj
xj+1
 L’integrale tra xj-1 e xj+1 della f(x) è dato da:
x j 1

x j 1
 f
x j 1
f ( x)dx 
2
j

 a x  x j   bx  x j  dx 
x j 1
x h
j
a
b
2 3

2
3




 2hf j  bh
 f j x  2 x  x j  3 x  x j 
3
x j h
x
Regola di Simpson (3)
 Calcoliamo quindi il coefficiente b:
 f ( x j  h)  f j 1  f j  ah  bh 2  f j 1  ah  bh 2  f j 1  f j



2
 f ( x j  h)  f j 1  f j  ah  bh  f j 1  ah  bh 2  f j 1  f j
 Sommando membro a membro le due equazioni si ha:
2bh  f j 1  2 f j  f j 1  b 
2
f j 1  2 f j  f j 1
2h 2
 Sostituendo il valore di b nella formula dell’integrale:
x j 1

2 3
4
1
1

5 ( 4)
f
(
x
)
dx

2
hf

bh

h
f

f

f

O
h
f
j
j

1
j
j

1
x


3
3
3
3

j 1

 Il valore dell’integrale calcolato con la formula di Simpson differisce dal
valore vero per un termine che è dell’ordine di h5 per la derivata quarta
della funzione calcolata in un punto (non noto) dell’intervallo [xj,xj+1]
 La formula di Simpson è esatta per polinomi fino al terzo grado
Formula di Simpson estesa
2h
a=x1
x2
2h
x3
x4
2h
x5
...
xN-2
xN-1
xN=b
x
 Dividiamo l’intervallo [a,b] negli (N-1)/2 intervalli:
 [x1,x3], [x3,x5], ..., [xN-2,xN]
 a=x1; b=xN
 h=(b-a)/(N-1)
 Si ha:
b

a
4
1  1
4
1 
4
1 
1
1
f ( x)dx  h f1  f 2  f 3   h f 3  f 4  f 5     h f N  2  f N 1  f N 
3
3  3
3
3 
3
3 
3
3
5
N  1  b  a  ( 4) 
4
2
4
2
4
1 
1
 1 

O 
 f   h  f1  f 2  f 3  f 4    f N  2  f N 1  f N   O 4 
2
3
3
3
3
3
3 
3
N 
 N  1 

Applicazione della regola di Simpson
 L’algoritmo è molto simile a quello usato per la regola del trapezio
 Anche in questo caso si procede per approssimazioni successive:
 nella n-esima iterazione vengono aggiunti 2n-2 valori di f(x) in
corrispondenza dei punti medi degli intervalli elementari definiti
dall’iterazione precedente
 alla n-esima iterazione l’intervallo di integrazione risulta diviso in 2n-1
intervalli elementari
 nella n-esima iterazione si sfrutta il risultato ottenuto con la regola del
trapezio nell’iterazione precedente
1
iterazioni
2
3
4
...
Implementazione dell’algoritmo (1)
 Poniamo ancora b-a=Δ
 Per la n-esima iterazione poniamo:
 Sn = integrale calcolato con la regola di Simpson
 Tn = integrale calcolato con la regola del trapezio
 Alla prima iterazione si ha:
S1  0
1
1

T1    f a   f b 
2
2

 notare che S1=0 perché nella prima iterazione si
considerano solo due punti, mentre per applicare la
regola di Simpson ne occorrono tre
Implementazione dell’algoritmo (2)
 Alla seconda iterazione si ha:
S2 
 1
4


f
a

2  3
3
 1
 1
4 


f  a    f b   T1 
f a  
2 3
23 
2

 3
 1
 1
 1
 


T2   f a   f  a    f b   T1  f  a  
2 2
2 2
2 
2

 2
 Alla terza iterazione si ha:
S3 

 1
4 
 2 
 4 
3  1




f
a

f
a


f
a


f
a


f
b







4  3
3 
4 3 
2 3 
4  3

1
 4
T2  
3
4 3
T3 
 4

f a   
4 3

 1
f a  

4 2
3 

f a 

4




f a   
4

1
 

T2   f  a   
2
4 
4


f a   
2

3 

f a 

4 


3  1

f a 
  f b  
4  2


Implementazione dell’algoritmo (3)
 Alla n-esima iterazione si avrà:
n2
1
 2 4
S n  Tn 1  n 1 
3
2 i 1 3
n2
1
 2
Tn  Tn 1  n 1 
2
2 i 1

 

f  a  n 1  i  1 n  2 
2
2 


 

f  a  n 1  i  1 n  2 
2
2 

 Poiché nel passare da un’iterazione alla successiva il
numero di punti utilizzati per il calcolo raddoppia, la
precisione del calcolo migliora di un fattore 16.
 Come nel caso precedente, la procedura iterativa viene
fermata quando:
Sn  Sn1   Sn1
dove ε è la precisione richiesta per il calcolo
Esempi (1)
 Supponiamo di voler calcolare l’integrale
1

dx

 0.785398
2
1 x
4
0
1
I 
con una precisione ε=10-6. Si ha:
Iterazione (n)
Tn
Sn
1
0.75
0
2
0.775
0.783333
3
0.782794
0.785392
4
0.784747
0.785398
5
0.785235
0.785398
Esempi (2)
 Supponiamo di voler calcolare l’integrale
2
1
I   dx  ln 2  0.693147 
x
1
con una precisione ε=10-6. Si ha:
Iterazione (n)
Tn
Sn
1
0.75
0
2
0.708333
0.694444
3
0.697024
0.693254
4
0.694122
0.693155
5
0.693391
0.693148
6
0.693208
0.693147
Calcolo di integrali impropri
 Casi di integrali impropri estesi ad intervalli finiti:
 la funzione integranda ha limite finito in uno degli estremi di
integrazione, ma non può essere calcolata (esempio sinx/x per x0)
 la funzione ha limite superiore + e/o limite inferiore -
 la funzione ha una singolarità integrabile in uno dei due estremi di
integrazione (esempio x-1/2 per x 0)
 la funzione ha una singolarità integrabile in un punto noto
dell’intervallo di integrazione (ci si può ricondurre al caso precedente)
 la funzione ha una singolarità integrabile in un punto non noto
dell’intervallo di integrazione (questo caso non verrà studiato)
 In tutti i casi in esame, per poter effettuare il calcolo, è necessario
che l’integrale esista e sia finito
 se l’integrale non esiste oppure è infinito, qualunque procedura di
calcolo sarà priva di senso e darà risultati errati!
 Studieremo il calcolo di integrali impropri con una singolarità in
corrispondenza di uno degli estremi di integrazione (o di entrambi)
 in tal caso non è possibile usare una formula chiusa, perché
implicherebbe il calcolo del valore della funzione integranda in
corrispondenza della singolarità
Regola del punto medio
 La regola del punto medio consiste nell’approssimare l’integrale
tra xj e xj+1 nel modo seguente:
x j 1
 f ( x)dx  hf
j 1 2
xj
ove f j+1/2 è il valore della funzione nel punto medio dell’intervallo
[xj,xj+1]
 se f(x)>0, l’integrale viene approssimato con l’area del rettangolo
di base h=xj+1-xj e altezza f j+1/2
y
fj
fj+1/2
fj+1
xj
xj+1/2
xj+1
x
Regola estesa del punto medio
h
a=x1
h
x3/2
x2
x5/2
h
x3
...
xN-1
xN-1/2 xN=b
 Dividiamo l’intervallo [a,b] negli N-1 intervalli:
 [x1,x2], [x2,x3], ... , [xN-1,xN]
 a=x1; b=xN
 h=(b-a)/(N-1)
 Si ha:
 f ( x)dx  h f
b
a
3/ 2
 f 5 / 2    f ( N 1) / 2

x
Applicazione della regola del punto medio (1)
 L’algoritmo è simile a quelli visti per la regola del trapezio
e per la regola di Simpson
 In questo caso, però, per poter usare il risultato ottenuto
dopo ogni iterazione come punto di partenza per
l’iterazione successiva, in ciascuna iterazione occorrerà
suddividere gli intervalli di partenza in 3 parti invece che
in 2
 questa complicazione nasce dal fatto che si utilizzano i
valori della funzione nei punti medi di ciascun intervallo
invece che negli estremi
Applicazione della regola del punto medio (2)
iterazioni
1
2
3
...
 nella n-esima iterazione vengono aggiunti 23n-2 valori di f(x)
in corrispondenza dei punti medi dei nuovi intervalli
 alla n-esima iterazione l’intervallo di integrazione risulta diviso
in 3n-1 intervalli elementari
 in ciascuna iterazione si sfrutta il risultato ottenuto
nell’iterazione precedente
Implementazione dell’algoritmo (1)
 L’integrale calcolato nella n-esima iterazione è:
n1
 3
I n  n 1 
3 k 1

 
1 
f a  n 1  k  
3 
2 

dove Δ=b-a
 Occorre cercare una formula ricorsiva che permetta di
legare In+1 a In
 Applicando la definizione precedente si ha:
I n 1 

3n
3n

k 1


f a  n
3

1 

k



2 

Implementazione dell’algoritmo (2)
 La sommatoria con l’indice k che varia da 1 a 3n può
essere spezzata in 3 sommatorie distinte:
 nella prima raggruppiamo i termini con k=1,4,7,...,3n-2
 in tali termini k=3k1-2 con k1=1,2,...,3n-1
 nella seconda raggruppiamo i termini con k=2,5,8,...,3n-1
 in tali termini k=3k2-1 con k2=1,2,...,3n-1
 nella terza raggruppiamo i termini con k=3,6,9,...,3n
 in tali termini k=3k3 con k3=1,2,...,3n-1



 3

5  3

3  3

1  
I n 1  n  f a  n  3k1     f a  n  3k 2     f a  n  3k3   
3 k1 1  3 
2  k2 1  3 
2  k3 1  3 
2  
3n1
3n1
3n1





 
5 
 
1 
 
1  
 n  f a  n 1  k1     f a  n 1  k 2     f a  n 1  k3   
3 k1 1  3 
6  k2 1  3 
2  k3 1  3 
6  
n1
1

 In  n
3
3
n1
 
 
5 
f
a

k

 
 

n 1 
3 
6 
k1 1  
3n1

 
1  
f a  n 1  k   
6  
 3 
n1
Implementazione dell’algoritmo (3)
 Nel passare da un’iterazione alla successiva il numero di
punti utilizzati per il calcolo viene triplicato, e dunque la
precisione del calcolo migliora di un fattore 9.
 Come negli altri casi, la procedura iterativa viene
fermata quando:
I n  I n1   I n1
dove ε è la precisione richiesta per il calcolo
 L’algoritmo sviluppato per il calcolo di integrali impropri
può essere usato anche per il calcolo di integrali
“normali”
 Il tempo di calcolo richiesto per un integrale improprio è
in genere maggiore rispetto a quello richiesto per un
integrale “normale” per via della singolarità
Esempio
 Supponiamo di voler calcolare l’integrale improprio
1
I   ln xdx  1
0
con una precisione ε=10-6. Si ha:
Iterazione (n)
In
1
-0.693147
2
-0.896696
3
-0.965451
4
-0.988471
5
-0.996156
6
-0.998718
7
-0.999573
8
-0.999858
9
-0.999953
10
-0.999984
11
-0.999995
12
-0.999998
13
-0.999999
14
-1.000000