6 oktober 2018, Brussel. Numerieke integratie en CAS. Guido Herweyers ex-docent KU-Leuven en VIVES Campus Brugge

1

21

^ste

T3 Vlaanderen Symposium, 6 oktober 2018, Brussel

Numerieke integratie en CAS Guido Herweyers

ex-docent KU-Leuven en VIVES Campus Brugge

1. Inleiding

Numerieke integratie heeft als doel een benaderende numerieke waarde te berekenen voor de bepaalde integraal ( )

b

a

f x dx



. Aan bod komen de rechthoekregels als bijzondere Riemann- sommen, de trapeziumregel en de regel van Simpson. We bespreken de verbanden tussen die regels en de foutafschattingen. Grafische voorstellingen en computeralgebra (CAS) bieden hier een meerwaarde, we illustreren dit met de TI-Nspire CX CAS.

2. De Riemann-integraal

Zij f een reële functie gedefinieerd op het interval [a,b]. Verdeel het interval [a,b] in n deelintervallen met lengte x_k x_k x_k1 , k 1, 2, : ,n

Kies in elk deelinterval



x_k1,x_k



een punt c_k (k 1, ,n) en bereken de bijhorende Riemann-som

1

( )

n

n k k

k

S f c x







 ^:

2 Laat n nu onbeperkt toenemen, waarbij de lengte van elk interval x_k (k1, 2, ) naar nul gaat, ,n dan convergeert S_n naar de Riemann-integraal ( )

b

a

f x dx 



 , als de functie f integreerbaar is op [a,b]. Een (stuksgewijs) continue functie op [a,b] is alvast integreerbaar.

Meetkundige interpretatie: ( )

b

a

f x dx



= “oppervlakte boven x-as - oppervlakte onder x –as”

3. Rechthoekregels

Verdeel het interval [a,b] in n deelintervallen met gelijke lengte b a

h x

n

    , dan geldt

xk   a k h , stel y_k  f x( )_k (k 0,1, ). ,n

a) Kies voor elk deelinterval



x_k1,x_k



het startpunt c_k x_k1 a



k1



h (k1, 2, ) , ,n dan wordt de Riemann-som de linkersom

 

0 1 1

0 1 1

1

( )

( ( 1) )

n

n n

k

l n y h y h y h

h y y y

h f a k h







      

   

 



  





linkersom met n=4

3 b) Kies voor elk deelinterval



x_k1,x_k



het eindpunt c_k x_k   a k h (k1, 2, ) , dan wordt de ,n Riemann-som de rechtersom

 

1 2

1 2

1

( )

( )

n

n n

k

r n y h y h y h

h y y y

h f a k h



      

   

 



 





rechtersom met n=4

c) Kies voor elk deelinterval



x_k1,x_k



het middelpunt ¹



¹2



2

k k

k

x x

c ^  a k h

    (k1, 2, ) , ,n

dan wordt de Riemann-som de middensom ¹₂

1

( ) ( ( ) )

n

k

m n h f a k h



 



  

middensom met n = 4

4 Observeer de convergentie van de rechthoeksommen naar ( )

b

a

f x dx



door n te laten toenemen:

4. De trapeziumregel

In plaats van te werken met oppervlakten van rechthoeken werken we nu met oppervlakten van trapezia (interpretatie als ( )f x 0voor ^x



^{a b,}



).Verdeel het interval [a,b] in n deelintervallen met gelijke lengte b a

h x

n

    , met x_k   a k h en y_k  f x( )_k (k 0,1, ). ,n

Door de grafiek van f op elk deelinterval te benaderen door een lijnstuk (lineaire interpolatie) door de punten (x_k1,y_k1) en (x y_k, _k) (k1, 2, ) , verkrijgen we de trapeziumregel: ,n

0 1 1 2 1

0 1 2 1

1

1

( ) ( ) ( ) ( )

2 2 2

( 2( ))

2

( ) ( ) 2 ( )

2

n n

n n

n

k

h h h

t n y y y y y y

h y y y y y

h f a f b f a k h









     

    

 

      











Trapeziumregel met n=4

5 5. De regel van Simpson

Verdeel het interval [a,b] in een even aantal deelintervallen n met gelijke lengte b a

h x

n

    ,

met x_k   a k h en y_k  f x( )_k (k 0,1, ). Benader de grafiek van f op elk opeenvolgend paar ,n deelintervallen door een parabool (kwadratische interpolatie) door de punten (x y_k, _k), (x_k1,y_k1) en (x_k2,y_k2) (k0, 2, 4,,n2)

Computeralgebra levert de oppervlakte onder de eerste parabool door de punten

0 1 2

( ,a y ) , (ah y, ), (a2 ,h y ) . Het eenvoudige resultaat is ( ₀ 4 _{1 2}) 3

h y  y y :

6 Sommatie van de oppervlakten onder alle parabolen levert de regel van Simpson:

0 1 2 2 3 4 4 5 6 2 1

0 1 3 5 1 2 4 6 2

/ 2 / 2 1

1 1

( ) ( 4 ) ( 4 ) ( 4 ) ( 4 )

3 3 3 3

( 4( ) 2( ))

3

( ) ( ) 4 ( (2 -1) )+2 ( 2 )

3

n n n

n n n

n n

k k

h h h h

s n y y y y y y y y y y y y

h y y y y y y y y y y

h f a f b f a k h f a k h

 

 



 

            

           

 

        



 





 

We definiëren nu eerst de verschillende numerieke integratiemethoden in een rekenmachine pagina van TI-Nspire (dit gaat vlot door kopiëren):

7 Vervolgens definiëren we f(x), a en b in een dynamische notitiepagina:

De verschillende methoden convergeren naar de exacte integraal (bereken lim

n).

Uiteraard levert s(n) de exacte integraal voor elke n bij integratie van een tweedegraadsfunctie.

Ook voor een derdegraadsfunctie levert de regel van Simpson steeds een exact resultaat (zie verder bij de studie van de fout)!

6. Verbanden tussen de numerieke integratiemethoden

Er zijn veel verbanden tussen de verschillende integratiemethoden, o.a. : ( ) ( )

( ) 2

l n r n

t n 

 (eenvoudig te bewijzen) 4 (2 ) ( )

( ) 3

t n t n

s n 

 ( ) ( )

(2 ) 2

t n m n

t n 

 (2n-strip precisie met n-strip berekening) ( ) 2 ( )

( ) 3

t n m n

s n 



8 We verifiëren deze eigenschappen met computeralgebra voor een concrete f(x), a en b:

Vergelijk ( ), ( ), ( )l n r n t n en ( )m n voor dalende en stijgende functies op een interval [a,b]:

Dalende functie, holle kant naar beneden:

bovensom ( )l n , ondersom ( )r n , ( ) ( )

( ) 2

l n r n

t n 

 ( )r n t n( )m n( )l n( ) Dalende functie, holle kant naar boven:

bovensom ( )l n , ondersom ( )r n , ( ) ( )

( ) 2

l n r n

t n 

 r n( )m n( )t n( )l n( )

9 Stijgende functie, holle kant naar beneden:

bovensom ( )r n , ondersom ( )l n , ( ) ( )

( ) 2

l n r n

t n 

 ( )l n t n( )m n( )r n( ) Stijgende functie, holle kant naar boven:

bovensom ( )r n , ondersom ( )l n , ( ) ( )

( ) 2

l n r n

t n 

 ( )l n m n( )t n( )r n( )

Algemeen geldt:

Voor een dalende functie op een interval [a,b] is de bovensom ( )l n en de ondersom ( )r n . Voor een stijgende functie op een interval [a,b] is de bovensom ( )r n en de ondersom ( )l n . Voor een functie met de holle kant naar boven op een interval [a,b] geldt ( )m n t n( ). Voor een functie met de holle kant naar beneden op een interval [a,b] geldt ( )t n m n( ).

10 Een schuifbalk voor n in de notitiepagina zorgt voor een dynamische aanpassing van de waarden:

Voor de trapeziumregel geldt ( ) ( ) ( )

b

t a

f x dxt n E n



, waarbij E n_t( )de fout is.

Voor een stijgende of dalende functie f op het interval [a,b] geldt steeds ( ) ( )

( ) 2

t

r n l n

E n 

 ,

aangezien ( ) ( )

( ) 2

l n r n

t n 

 enomdat de exacte waarde van de integraal gelegen is tussen de ondersom en de bovensom (de waarden ( )l n en ( )r n ).

Nu geldt r n( )l n( )h y



1y2y_n



h y



0y1y_n1



h y( _n y0) h ( ( )f b  f a( )) ,

zodat ( ) ( ) ( )

t 2

E n h f b  f a of ( ) ( ) ( )

t 2

E n b a f b f a n

   .

11

Meetkundige voorstelling van ( )r n l n( )voor een stijgende functie op [a,b]

7. Fouten ( )E n bij numerieke integratiemethoden

a) Linkersom of rechtersom

Als 'f continu is op het interval [ , ]a b en f x'( ) M voor elke x[ , ]a b dan geldt

( ) ( ) ( )

b

r a

f x dxr n E n



^met

( )2

( ) 2

r

E n M b a n

  , ( ) ( ) ( )

b

l a

f x dxl n E n



^met

( )2

( ) 2

l

E n M b a n

 

b) Middensom

Als ''f continu is op het interval [ , ]a b en f ''( )x Mvoor elke x[ , ]a b dan geldt

( ) ( ) ( )

b

m a

f x dxm n E n



^met

3 2

( )

( ) 24

m

E n M b a n

  (omgekeerd evenredig met n²)

c) Trapeziumregel

Als ''f continu is op het interval [ , ]a b en f ''( )x Mvoor elke x[ , ]a b dan geldt

( ) ( ) ( )

b

t a

f x dxt n E n



^met

3 2

( )

( ) 12

t

E n M b a n

  (omgekeerd evenredig met n²)

d) Regel van Simpson

Als f⁽⁴⁾ continu is op het interval [ , ]a b en f^{( 4)}( )x M voor elke x[ , ]a b dan geldt

( ) ( ) ( )

b

s a

f x dxs n E n



^met

5 4

( )

( ) 180

s

E n M b a n

  (omgekeerd evenredig met n⁴)

12 Merk op dat de regel van Simpson de exacte integraalwaarde levert voor alle veeltermen tot en met graad 3, maar niet voor elke vierdegraadsveelterm (illustreer dit met een voorbeeld); de regel van Simpson heeft nauwkeurigheidsgraad 3.

Als voorbeeld schatten we het minimum aantal deelintervallen n die nodig zijn om

2

1

1dx ln(2) 0.69314718055989

x  



te benaderen met de regel van Simpson zodat E n_s( ) 10^6:

bepaal n zodat

5 6 4

( )

( ) 10

s 180 E n M b a

n

 

  . Voor f x( ) ¹ x ¹

x

   vinden we f^{( 4)}( )x 4!x ^{5 24}₅ x

   .

Een bovengrens van f⁽⁴⁾( )x op het interval [1, 2] is M 24.

Uit

5 6 4

( )

180 10 M b a

n

 

 of ²⁴₄ 10 ⁶ 180n

  volgt ^{4 2} ₆ 15 10

n  _

 zodat ⁴ 2 ₆

19,11 15 10

n _ 

 .

Als schatting vinden we het even getal n 20. Deze schatting is pessimistisch want de correcte minimale waarde is n 14.

Voor de trapeziumregel vinden we voor n 14:

( ) ( ) ( )

t 2

E n b a f b f a n

   (p. 10 voor een dalende functie) zodat 2 1 1 1

(14) 1 0, 018

2 14 2 56

Et 

    

 .

Anderzijds geldt

3 2

( )

( ) 12

t

E n M b a n

  met f ''( )x Mvoor elke x [1, 2] _(p.11). Nu is f ''( )x ²₃

 x , met M 2 vinden we een betere afschatting

3

2 2

(2 1) 6

(14) 2 0, 00085

12 14 14

Et 

  

 ^.