Det bestemte integralet med rektangelmetoden

Vi skal finna ein tilnærma verdi for det bestemte integralet mellom 1 og 7 av denne funksjonen:

$f(x) = x³ - 9x² + 23x - 5$.

Sidan f(x) > 0 i heile dette området, er det det same som å finna arealet som er markert under. Hvis du har hatt om bestemte integraler, så kan du rekna ut dette sjøl. Fasiten er 96.

Det enklaste er å bruka rektangelmetoden. Det betyr at vi vil finna ein tilnærmingsverdi for integralet ved å summera eit antal rektangel slik figuren under viser. Her har vi delt opp i N = 12 rektangel:

La oss begynna med å finna arealet av rektangel nr. 5, som vist i grønt på figuren. Vi legg merke til at høyden av rektangelet er lik funksjonsverdien på venstre side, dvs. f(3), så vi definerer først funksjonen vår og finn denne verdien:

def f(x):
    return x**3 - 9*x**2 + 23*x - 5

print("f(3)=",f(3))

f(3)= 10

Så kan vi leggja inn nedre grense a ( = 1) og øvre grense b ( = 7) for arealet, samt antalet rektangel N ( = 12), Dermed kan vi rekna ut bredden av rektangela, som vi kallar dx:

a = 1 # nedre grense
b = 7 # øvre grense
N = 12 # antal rektangel
dx = (b-a)/N # Bredden av rektangela

print("bredden av rektangela, dx =",dx)

bredden av rektangela, dx = 0.5

Det vi no treng å gjera er å laga ei løkke som reknar ut arealet til kvar av rektangela, og summerer disse. Arealet av eit rektangel er gitt ved $f(x)\cdot dx$. Men vi treng først ein generell formel som vi kan bruka inne i løkka, nemlig x = a + n*dx Her er n det som skal bli tellevariabelen i løkka. Og før vi lagar løkka ,vil vi testa om denne formelen fungerer for n = 4. Husk at vi startar å telja med 0, slik at dette er det femte rektangelet.

n = 4 # Dette legg vi inn for testformål. Denne linja skal ikkje med i den endelige versjonen!
x = a + n*dx # x-verdien er lik nedre grense (1) + 4 ganger 0,5
h = f(x)
A = dx*h
print("x =",x,"gir høyden",h,"og arealet:",A)

x = 3.0 gir høyden 10.0 og arealet: 5.0

Det ser ut til å fungera, for 10 * 0.5 er jo som kjent lik 5! Så alt vi treng å gjera no er å legga denne koden inn i ei løkke, og summera for kvar iterasjon. Sumvariabelen vår er S, og vi lar no tellevariabelen n gå frå 0 til N:

S = 0 # Før løkka lar vi S få startverdien 0. 
for n in range(N):
    x = a + n*dx
    h = f(x)
    A = dx*h
    S += A # Her legg vi til A for kvar iterasjon.

print("Summen av",N,"rektangel er:",S)

Summen av 12 rektangel er: 84.75

For å sjekka om arealet nærmar seg det riktige svaret, lagar vi ein ein funksjon som vi kallar rektangelmetoden(f,a,b,N), som vi kan kjøra med ulike N-verdiar. Vi testar den først ned N = 12, for å sjå om den fungerer:

# Funksjon som reknar ut summen av N rektangel frå a til b for funksjonen f:
def rektangelmetoden(f,a,b,N):
    S = 0
    dx = (b-a)/N #Bredden av rektangela
    for n in range(N):
        x = a + n*dx
        h = f(x)
        A = dx*h
        S += A
        
    print("Summen av",N,"rektangel mellom a =",a,"og b =",b,"er:",S)

# Så kjører vi funksjonen vår med verdiane under:
rektangelmetoden(f,1,7,12)

Summen av 12 rektangel mellom a = 1 og b = 7 er: 84.75

La oss no testa rektangelmetoden med høgare N-verdiar:

rektangelmetoden(f,1,7,12)
rektangelmetoden(f,1,7,100)
rektangelmetoden(f,1,7,1000)
rektangelmetoden(f,1,7,10000)

Summen av 12 rektangel mellom a = 1 og b = 7 er: 84.75
Summen av 100 rektangel mellom a = 1 og b = 7 er: 94.5708
Summen av 1000 rektangel mellom a = 1 og b = 7 er: 95.85610800000002
Summen av 10000 rektangel mellom a = 1 og b = 7 er: 95.9856010800001

No ser det ut til at vi har eit fungerande program. Men vi vil gjerne også at det tegnar funksjonen og alle rektangela. Samtidig legg vi inn ein variabel offset som blir forklart etterpå. Men først heile koden:

import matplotlib.pyplot as plt
import numpy as np

def f(x):
    return x**3 - 9*x**2 + 23*x - 5

def rektangelmetoden(f,a,b,N,offset = 0): # Hvis vi ikkje gir offset, så bruker vi verdien 0
    S = 0
    dx = (b-a)/N #Bredden av rektangela
    
    plt.xlabel('x')
    plt.ylabel('y')
    plt.title("Rektangelmetoden:")
    plt.plot([a,b],[0,0],color='blue') # Horisontal linje frå a til b

    for n in range(N):
        x = a + n*dx
        h = f(x +offset*dx)
        A = dx*h
        S += A
        
        plt.plot([x,x],[0,h],color='blue') # Vertikal linje frå (x,0) til (x,h), dvs venstre side.
        plt.plot([x+dx,x+dx],[0,h],color='blue') # Vertikal linje i (x+dx,0) til (x+dx,h), dvs. høgre side.
        plt.plot([x,x+dx],[h,h],color='blue') # Horisontal linje frå (x,h) til (x+dx,h)

    x = np.linspace(a,b,100) # Genererer 100 x-verdiar mellom a og b       
    plt.plot(x,f(x),color='red') # Plotter sjølve funksjonen
    
    plt.show()
    
    print("Summen av",N,"rektangel mellom a =",a,"og b =",b,"er:",S)

rektangelmetoden(f,1,7,12)

Summen av 12 rektangel mellom a = 1 og b = 7 er: 84.75

Som figuren viser, er høyden til alle rektangela rekna ut frå x-verdien på venstre side. Men variabelen offset gir oss mulighet til å rekna den ut eit lite stykke til høgre, som er gitt ved offset*dx. Som default er verdien 0. Men la oss testa kva vi får når vi set offset lik 1:

rektangelmetoden(f,1,7,12,1)

Summen av 12 rektangel mellom a = 1 og b = 7 er: 108.75

Vi ser at no blir høyden rekna ut frå x-verdien på høgre side. Men dette gav oss ikkje betre svar. Vanligvis vil ein offset-verdi på 0.5 fungera bra. Det vil sei at vi reknar høyden som funksjonsverdien av midtpunktet i rektangela:

rektangelmetoden(f,1,7,12,0.5)

Summen av 12 rektangel mellom a = 1 og b = 7 er: 95.625

No er vi mykje nærare det rette svaret. Som siste test, kan vi la N vera 100:

rektangelmetoden(f,1,7,100,0.5)

Summen av 100 rektangel mellom a = 1 og b = 7 er: 95.99459999999999

Det er alltids meir vi kan leggja inn i koden, som feks innlesing av verdiane a, b, N og offset, eller kontroll av disse verdiane mm. Men dette overlet eg til lesaren!

Oppgaver

1. Prøv programmet med andre verdiar for N.
2. Prøv med andre grenser a og b. Rekn ut det rette svaret, eller bruk GeoGebra-kommandoen "Integral(Funksjon,Start,Slutt)"
3. Sjekk også med andre funksjonar. Finn metoden vår gode svar?