TECHNISCHE HOGESCHOOL DELFT
AFDEUNG DER SCHEEPSBOUW- EN SCHEEPVAARTKUNDE CENTRALE WERKGROEP WISKUNDERapport No. 3.
BEPALING VAN HET OPPERVLAK VAN EEN WILLEKEURIGE
FUNKTIE
y =
f(x)
DOOR RECURSIEF GEBRUIK VAN
EEN INTEGRATIE PROCEDURE OP VERSCHILLENDE
DEELINTERVALLEN.
Ing. A.P.de Zwaan
augustus 1977
DeRI Unversi1y of Technology Ship Hydromechanics Laboratory Mekelweg 2
Delft 2208
INHOIJD
Inleiding.
1.2 Formules.
Beschrijving van het programma. 2.1 Algemeen
2.1.1. Programma gegevens. 2.1 .2. Motivering.
2.1.3. Opzet.
2.2 Organisatie van het prograliuiia. 2.2.1. Gebruikte methoden. 2.2.2. Formules.
2.2.3. Benodigde randapparatuur. 2.2.4. Verkiaring der variabelen. 2.3 Bestandsorganisatie.
2.3.1. Gebruik van de.procedure. 2.3.2. Gebruiksmogelijkheden. 2.3.3. Invoer.
2.3.4. Blokdiagram.
2.3.5. Listing van de procedure.
blz.: 2 2 6 6 6 6 6 6 6 6 6 6 7 7 7 7 8
1. INLEIDING.
De publikatie geeft aan, hoe van een "willekeurige" funktie y
= f(x),
tussenbepaalde grenzen en met een opgegeven nauwkeurigheid, het opperviak berekend wordt. De methode berust op de regel van Simpson met de korrektieterm van Richardson, De opzet is zodanig, dat stukken van de funktie die eeu te grote fout geven
worden verfijnd en de rest van de funktie die geen te grote fout geeft, niet wordt
verfijnd.
Tevens wordt de nauwkeurigheid zodanig bijgestuurd, dat als een bepaald gedeelte van de funktie met een grotere nauwkeurigheid wordt berekend dan opgegeven, dit
in mindering wordt gebracht voor het volgende te berekenen gedeelte van de
funktie.
Dit heeft bet voordeel, dat hele stijle gedeelten in de funktie met een kleinere nauwkeurigheid magen worden berekend, als men er voor zorgdraagt dat die stukken het eerst worden berekend, die het dichtst bij de opgegeven nauwkeurigheid liggen.
Zou deze aanpassing niet toegepast worden, dan bestaat de kans, dat de fout in het het oppervlak van de hele stijle gedeelten van de funktie, ondanks de grootste verfijning, groter blijft dan de opgegeven nauwkeurigheid.
1.2. Formules,
De hoofdformule van Simpson is;
B (f(0) + 4 f(2) f(4)) (B-A)5 (4) B-A S1 = J f(x)dx 6 2880 f(t1) A
met t1 uit[A,B zio fig. i
B
3
Figuur i..
Door verkleining van het interval door de punten i en 3 tussen te voegen is:
B (AB)/2 B S2 =
f f(x)dx
=J
f(x)dx +f f(x)dx
= A A (A+B)/2 (4) (4) B-A (f( 0) + 4f(1) + 2f(2)+ 4f(
3)+ f(4))
((B-A)/2)5f(t)
+ f(t2) 1440 2M-Bi
.iA+B
met t1 uit [A,
2
J
eri t7uit
L' 'J
(f(0)
+ 4f(1) + 2f(2) + 4f(3)+ f(4))
((B-A)/2)5 (4) =f
f(X)dx
B-A f (t) 12 1440 Amet t uit [A,BJ
Stel B-A 4 B S1 =
J f(x)dx
A B 1 S2J
f(x)dx A f(4)) Ï6 h5 (4) f(t1)
15(4)
2f(2) + 4f(3) + f(4)) - h f(t)We kunnen zeggen S1 = S2 dus S1 - S2 = O (met korrektieterm). Hieruit voigt: 16 5 (4) (4) 1/3 h (f(0) - 4f(1) + 6f(2) - 4f(3) + f(4)) - h r
(t1)
h5 f(t)
= O (4) (4) Stel f (tj) f (t): 1/3 h (f(0) - 4f(1) +6f(2) - 4f(3) + f(4)) 1/3 h5 f4(t) ett
uit ,B] Hieruit voigt: de korrektie voor S2 iS h (f(0) - 4f(1) 6f(2) - 4f(3) + f(4))de korrektie voor S1 is h (f(0) - 4f(i) + 6f(2) - 4f(3) + f(4))
De korrektie voor S1 is dus 16 maal zo groot.
Bepaling van de korrektieterm van Richardson.
A Sj + p S2
Zodeen lineaire kombinatie
+ zodanig dat voor een n-de graads
f(x) exact geldt. p
B
J f(x)dx
+ p S2 met XA+p
A
De lineaire kombinatie voor de regel van Simpson is:
B
J f(x)dx
42S2 - Si
-1
A
Deze lineaire kombinatie tussen Si en S2 geeft een kleinere fout dan voor
elk van deze afzonderlijk.
De fout in S1 is 16 E en in S2 i c E. B J f(x)dx 16S7 - Si - S2 Sj - S2 - s h5 f
(t)
(4) IS 15 2 2 A of B f(x)dx = S2 - h (f(0) - 4f(1) + 6f(2) - 4f(3)+ f(4))
A -3-= 2/3 h (f(0) + 4f(2) = 1/3 h ( f(0) 4f(1)-4--De korrektie voor Sj2 is
Kl = h (f(0) - 4f(1) + 6f(2) - 4f(3)
+ f(4))
Dit heet de "Korrektieterm van Richardson't.
Als geldt (f(0) - 4f(1) + 6f(2) - 4f(3) + of 180 I
f(0) -
4f(1) + 6f(2) - 4f(3) + <danis: I1S2 1<
sdus zeker ¡I - (Sj2
- -j--h- T) I s
met T 1/3 h jf(0) - 4f(1) 6f(2) - 4f(3)
+ f(4)J.
De gedachte, die hierachter zit is, dat de stukken die een te grote fout geven worden verfijnd en de rest die geen te grote fout geven niet behoeven te worden
verkleind. Toepassing. De exakte oplossing is B f(x)dx =
I
h B-A 4 [.2.3 ADe oplossing met de regel van Simpson inclusief de korrektieterin van Richardson is:
S 1/3 h (f(x) + 4f(xj+h) + 2f(xj+2h) + 4f(x1+3h) + f(x1+4h)) +
-
h (f(x) - 4f(x+h) + 6f(xí+2h) - 4f(xj+3h) + f(x+4h)) =0 - K
[i.2.]Geist wordt: I-S < s
Gesteld wordt dat bet interval [A,B] in n stukken moet worden opgedeeld voordat
aan deze eis is voldaan.
Voor elk stuk wordt het opperviak berekend volgens [1.2.4)
Het totale berekende opperviak moet voldoen aan:
EI-E
s Resumerend A I I =J
f(x)dx 0 1s1
2/3 h (f(0) + 4f(2)+ f(4))
S2
= 1/3 h (f(0) + 4f(1) + 2f(2) B I 2 3 4 + 4f(3) f(4))'b-'-'
-Het gehele rekenproces met aangepaste nauwkeurigheid verloopt als voigt:
E = n K-1
= --
C - E lJj voor -c = 1, 2, n [1.2.5] i= IVoigens de vergelijking [i.2.4] is: S O - Kn iri [Xj,X]
= h.(f(x1) - 4f(x1+h) + 6f(x2h) - 4f(xi+3h) + f(x)) -- h.T B-A met h 1 2 c+1 I met 'Sk+1 n + ('Skik) - Pk 4.n 180 Neem cS B-A C en T I B-A 180 dan is K 4n B-A C Kn
-Als dus voor opperviak Sk, voor k = 1, 2, 3 ...n, de eis wordt gesteid
180 . C
Tk = --- C, dan voigt daaruit, dat I -
sj
dus oak - < c.
Belangrijk is, dat de afschatting van Tk C onafhankelijk is van het aantal
deelintervallen.
Het voordeel hiervan is, dat men tijdens het rekenproces kan bepalen of een bepaaid
stuk van de funktie verfijnd moet worden of niet.
Men zit dus niet aan een bepaalde verdeling vooraf vast.
In het rekenproces wordt 'S voor opperviak Sj aangepast zoals aangegeven in verge-iijking [1.2.5]
-5--nStel de fout in SK is k dan moet gelden Z
pj
c.1=1
I1k-Ski 'Sk, neern 'S - en k 'Si dan wordt aan deze eis voldaan.
De fout in Sk is p, deze fout macht zijn 'Sk.
Het volgende opperviak
S+j
mag nu met een kleinere nauwkeurigheid worden berekendterwiji aan de eis + - (Sj + S+1) - c toch voldaan wordt.
2. BESCHRIJVING VAN ITET PROGRAMMA. 2.1. Algeineen. 2.1.1 Programma geg Taal: Geheugen: Rekentijd: Naam: 2.1.2 Motivering.
In de lijst van numerieke subprogramma's van het Rekncentrum is de
procedure "INSIG" opgenomen.
Deze procedure is goed bruikbaar voor funkties met een "normaal"
verloop.
Voor funkties, die een "heel stiji" verloop hebben is deze procedure niet bruikbaar (kan in veci gevallen de opgegeven nauwkeurigheid niet bereiken en gebruikt zodoende te veci geheugen).
Daarom is de procedure "IREC" geschreven.
Deze procedure kan voor elke kontinue nwaardige funktie het opper-viak bepalen met een opgegeven absolute nauwkeurigheid, ongeacht de
stij iheid.
2.1.3 Opzet.
Het programma is in de vorm van een recursieve procedure geschreven. De procedure bepaalt hoe en waar de funktie verfijud wordt.
De nauwkeurigheid wordt in het rekenproces steeds bijgestuurd, zodanig dat, de som van de opperviakken verminderd met de werkelijke som zo dicht mogeiijk ligt bij de opgegeven nauwkeurigheid.
2.2. Oranisatie van het proßramma.
2.2.1 Gebruikte methoden.
De methode voor dc bepaling van het opperviak berust op de methode van Simpson met de korrektie van Richardson.
2.2.2 Formules.
Zie hoofdstuk I
2.2.3 Benodigde randapparatuur. Geen.
2.2.4 Verkiaring van de variabelen. evens. Algol-60 variabel < I min Procedure "IREC".
-6-[ Eenheid Omschrijving Mathematisch Programma x X - De onafhankelijkverander-lijke van de funktie y f(x) f(x) Y - De funktie f(x) : De integrand
A A - De ondergrens van het
inter-gratie interval
B B - De bovengrens van het
inte-gratie interval
C EPS - De absolute nauwkeurigheid
waarmee de integraal herekend moet worden T T I-1 T2 Dc korrektie vanRichardson op de konstante 1/45 h na.
2.3. BestandsorBanisatie.
2.3.1 Gebruik van de procedure.
Aan de dekiaraties van het programma rnoet worden toegevoegd: 'REAL' 'PROCEDURE' IREC; 'CODE';
De procedure staat in de standaardbibliotheek SBAL.LIBCWW van de "Centrale Werkgroep Wiskunde" van de afdeling der Scheepsbouw- en Scheepvaartkunde.
2.3.2 Gebruiksmogelijkheden.
Enkele voorbeelden voor het gebruik
- tNT: = IREC (x, x2 + 2x + 1, 0,
10, 1O)
- A: = ; B: = ; EPS: =
INT: = IREC (x, x2 + 2x 1, A, B, EPS)
- tNT: = IREC (x, y(x), A, B, EPS) Hierin is y(x) een real procedure.
Opmerking: Indien men de foutmelding krijgt:
"Too many calls of procedures", dan is de waarde EPS te
klein genomen. 2.3.3 Invoer.
De kop van de procedure ziet er als voigt uit:
Real Procedure IREC (X,Y,A,B,EPS) Value A,B,EPS; Real X,Y,A,B,EPS; Invoer: Y = f(x) De integrand
A De ondergrens van bet integratie interval
B De bovengrens van het integratie interval EPS De absolute nauwkeurigheid waarmee de integraal
berekend moet worden
-Uitvoer: IREC Het gevraagde opperviak
Symbool
Eenheid Omschrijving
Mathematisch Programma
DELTA - T < DELTA
DELTANAX = = EPS
Si I - Deeloppervlak van de funktie
f(x) op de konstante 1 na 12
S IREC - Ret totale oppervlk:
0,1 2.3.4. Blokdiagram.
-8-xO = a x4 = b Bereken; T voor [xO,x4]L1
De waarde van X4 J' y dx wordt op de sommatie-ClijO
gep1aattI
Blad 2: Blokdiagram.
pl ,qI
(p2 ,q2)
betekent:
DELTA aanpassen voor bet volgen-de integratie interval p2,q2
De waarde van ydx wordt op
desommatielijn geplaatst.______
--I -9-t y Bereken: TI voor Lap1(b-a),a+q1(ba) T2 voor La+p2(b-a),aq2(ba)] xO a+pI(b-a) xO = a+p2(b-a) x4 a+ql(b-a) x4 = a+q2(b-a) T = TI T = T2L
DELTA aanpassen voor het
volgen-de integratie interval [pi,qiJ
X4
De waarde van J ydx wordt op xO
Blad 3: Blokdiagram.
io
-_______
I
o., 1/8(l//)4
o,'/
J
jz)
(O½4)
.-k
o H
(Y'a,/8
L----j 7/g(7/8,4)
F(o, ',/3')---
L---2.3.5
Listin
vn d
prnccHrtre
00010 'REAL' '
