Transcript Beräkningar i naturvetenskap och teknik
Numeriska beräkningar i Naturvetenskap och Teknik
Dagens ämne: Approximationsproblemet Minstakvadratmetoden Interpolationsproblemet Polynomanpassning Splines
Numeriska beräkningar i Naturvetenskap och Teknik
Ett exempel
T
20 30 40 50 60 70 80
L
f
(
T
) 48 .
56 50 .
19 50 .
60 50 .
71 51 .
76 53 .
48 54 .
83 54 52
L
50 48
Modell
L
f
(
T
)
c
0
c
1
T
20 40 60 80
T
Varför avviker mätvärdena från modellen om modellen är korrekt?
Numeriska beräkningar i Naturvetenskap och Teknik
Hur bestämma den ‘bästa’ linjen?
Modell
L
f
(
T
)
a
bT
54 52
L
50 48 20 40 60 80
T
Numeriska beräkningar i Naturvetenskap och Teknik
Avstånd mellan linje och mätpunkter...
54 52
L
50 48
d
1
d
2
d
3
d
4
d
5
d
6 20 40 60
d
7 80
T
Numeriska beräkningar i Naturvetenskap och Teknik
Norm Hur minimera avstånd mellan linje och mätpunkter?
d i
Största avvikelsen så liten som möjligt Approximation i maximum norm
d i
2
Summan av avvikelserna i kvadrat så liten som möjligt Approximation i Euklidisk norm Enklare att beräkna!
Numeriska beräkningar i Naturvetenskap och Teknik
Matrisformulering: Ett exempel
f x
3 1 0 2 3
f
(
x
) 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9
f
(
x
)
c
0
c
1
x c
0
c
0
c
0 3
c
1
c
1 2 0 .
9 .
1 0 .
6
c c
0 0 2
c
1 3
c
1 0 .
6 0 .
9
med
Ac
f
Överbestämt ekvationssystem!
1 1 1 3 1 0
c c
1 0 2 .
1 0 .
9 0 .
6 1 2 0 .
6 1 3 0 .
9
x
Numeriska beräkningar i Naturvetenskap och Teknik
Matrisformulering: Ett exempel
1 1 1 1 1 3 1 0 2 3
A T Ac
A T f c c
0 1 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9 1 3 1 1 1 0 1 2 1 3 1 1 1 1 1 3 1 0 2 3
c
0
c
1 1 3 1 1 1 0 1 2 1 3 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9
Numeriska beräkningar i Naturvetenskap och Teknik
Matrisformulering: ett exempel
A T Ac
A T f
5 1 1 23 c 0 c 1 2 .
1 11.1
c 0 c 1 0 .
521 0.505
f
(
x
) 0 .
521 0 .
505
x
Numeriska beräkningar i Naturvetenskap och Teknik
Matrisformulering: ett exempel
f
(
x
) 0 .
521 0 .
505
x
Numeriska beräkningar i Naturvetenskap och Teknik
Generell formulering: Beroende på modell kan mätdata förstås besrkrivas av andra funktionsuttryck än den räta linjen. I generella termer söks en funktion f* som approximerar f’s givna värden så bra som möjligt i euklidisk norm. Specifikt, ovan söktes en lösning
f
(
x
)
c
0
c
1
x
men man kunde ha sökt en lösning av annan form (ev. för andra data)
f
(
x
)
c
0 sin(
x
)
c
1 cos(
x
)
etc...
Numeriska beräkningar i Naturvetenskap och Teknik
Man kan således allmänt skriva:
f
(
x
)
c
0 0
c
1 1 ...
c n
n
f(x) är mao en linjärkombination av givna funktioner
0 , 1 ,...,
n
där koefficienterna
c
0 ,
c
1 ,...,
c n
söks. Man kan i likhet med ett vektorrum se det som att funktionerna
0 , 1 ,...,
n
spänner upp ett funktionsrum (ett rum av denna typ som uppfyller vissa villkor kallas ett Hilbert rum, jmf. kvantmek)
Numeriska beräkningar i Naturvetenskap och Teknik
I fallet med den räta linjen så är ex.vis
0 1 1
x
I en geometrisk jämförelse så spänner dessa två funktioner, som kan ses som två vektorer i funktionsrummet, ett plan U:
f
sökt funktion ”vektor” 0
0 1
f
*
approximationsfunktion ”vektor” 1
1
x
Minsta avståndet från planet ges av en normal . Den minsta avvikelsen mellan f* och f fås då f*-f är ortogonal mot planet U!
Numeriska beräkningar i Naturvetenskap och Teknik
Normalekvationerna Eftersom vi är intresserade av anpassning av m st mätvärden så kan vi lämna bilden av det kontinuerliga funktionsrummet och betrakta f(x) som en m dimensionell vektor med värdena:
f
f
(
x
1 )
f
(
x
2 ) ...
f
(
x m
)
Som skall uttryckas mha
0 0 ( 1 ) 0 ( 2 ) ...
0 (
m
)
samt
1 1 ( 1 ) 1 ( 2 ) ...
1 (
m
)
För den räta linjen:
0 1 1 ...
1
x
1 1
x
1 ...
x m
Numeriska beräkningar i Naturvetenskap och Teknik
Ortogonalitetskravet ger nu ekvationerna:
(
f
(
f
*
f
*
f
, 0 ) , 1 ) 0 0
med
f
*
c
0 0
c
1 1 (
c
0 0 (
c
0 0
c
1 1
c
1 1
f f
, 0 , 1 ) ) 0 0
c
0
c
0 ( 0 , 0 ( 0 , 1 ) )
vilket ger Normalekvationerna:
c
0 ( 0
c
0 ( 0 , , 0 ) 1 )
c
1 ( 1 1 , 0
c
1 ( 1 , 1 ) ) (
f
(
f
, 0 ) , 1 )
c
1 ( 1 1 , 0 )
c
1 ( 1 , 1 ) (
f
(
f
, 0 ) , 1 ) 0 0
Numeriska beräkningar i Naturvetenskap och Teknik 0 1
Normalekvationerna:
c
0
c
0 ( 0 ( 0 , 0 , 1 ) )
c
1 ( 1 , 0
c
1 ( 1 , 1 ) ) ( (
f f
, 0 , 1 ) ) ( 0 ( 0 , , 0 1 ) ( 1 ) ( 1 , , 0 1 ) )
c
0
c
1 (
f
(
f
, 0 ) , 1 ) .
.
.
.
.
.
0 1 .
.
.
.
.
.
c
0
c
1 0 1 .
.
.
.
.
.
f
Numeriska beräkningar i Naturvetenskap och Teknik
Tillbaka till exemplet: Data: Modell:
f
(
x
)
c
0
c
1
x x
3 1 0 2 3
f
(
x
) 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9 1 1 1 1 0 1 1 ( 3
x
) 1 0 2 3 A T 1 1 1 1 1 3 1 0 2 3 A 1 3 1 1 1 0 1 2 1 3 c
c
0
c
1 A T 1 1 1 1 1 3 1 0 2 3 f 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9
Numeriska beräkningar i Naturvetenskap och Teknik
Slutsats: Givet mätdata
f
under antagande av modellen
f
*
c
0 0
c
1 1
c
2 2
c
3 3 ...
c n
n
||
f
*
f
|| 2 0 , 1 , 2 , 3 ...
n
då
f
*
f
minimeras är ortogonal mot basvektorerna Koefficienterna c 1 , c 2 , c 3 , c n bestäms ur
Numeriska beräkningar i Naturvetenskap och Teknik
Normalekvationerna
c
0 ( 0
c
0 ( 0 , 0 , 1 ) )
c
1 ( 1 , 0
c
1 ( 1 , 1 ) ) ...
...
c n c n
( (
n n
, 0 ) , 1 ) ( (
f f
, 0 , 1 ) ) ...
c
0 ( 0 ,
n
)
c
1 ( 1 ,
n
) ...
c n
(
n
,
n
) (
f
,
n
)
eller
A T Ac
A T f
där kolumnerna i A ges av:
0 , 1 , 2 , 3 ...
n
Numeriska beräkningar i Naturvetenskap och Teknik
Notera 1: Funktionerna
0 , 1 , 2 , 3 ...
n
(jämför vektorer i ett vektorrum) måste vara linjärt oberoende Notera 2: Antag att vårt tidigare problem hade varit (x koord -996)
x
999 997 996 994 993
f
(
x
) 2 .
1 0 .
9 0 .
6 0 .
6 0 .
9 5 4979 4979 4958111 c 0 c 1 2 .
1 2102.7
jmf med
5 1 1 23 c c 0 1 2 .
1 11.1
Numeriska beräkningar i Naturvetenskap och Teknik
Gausseliminering i korthet:
x
1
x
1
x
1
x
1 2
x
2 4
x
2 8
x
2 16
x
2 3 9
x x
3 3 27
x
3 4
x
4 16
x
4 64
x
4 2 10 44 81
x
3 256
x
4 190 1 2 3 4 2 1 4 9 16 10 1 8 27 64 44 1 16 81 256 190
1 1 Numeriska beräkningar i Naturvetenskap och Teknik
Gausseliminering i korthet:
1 2 3 4 2 1 4 9 16 10 1 8 27 64 44 1 16 81 256 190 1 1 1 3 7 1 2 3 4 | 2 0 2 6 12 | 8 0 0 6 24 | 18 0 0 36 168 | 132 1 1 2 0 2 0 6 3 4 | 2 6 12 | 8 24 60 | 42 0 14 78 252 | 188 1 1 3 1 7 6 1 2 3 4 | 2 0 2 6 12 | 8 0 0 6 24 | 18 0 0 0 24 | 24
Numeriska beräkningar i Naturvetenskap och Teknik
Interpolationsproblemet Approximationen till mätdata antas gå genom datapunkterna, dvs man har tilltro till att felen i mätdata är små.
f
Linjär interpolation
(
x
)
f
0
f
1
x
1
f x
0 0 (
x
x
0 )
Alt vid ekvidistanta data
f
(
x
)
f
(
x
0
ph
)
f
0
p
f
0
x
f f 2 f 3 f 0 1.80
1.03
0.5
2.83
0.33
1.36
0.11
1.0
4.19
0.44
1.80
1.5
5.99
x f
0 0
ph h f
1
x
1
Numeriska beräkningar i Naturvetenskap och Teknik
Felkällor 1. Mätdata, E f
f E f x
0
x x
1
Numeriska beräkningar i Naturvetenskap och Teknik
Felkällor 2. Trunkationsfel Dessa vore noll för ett exakt förstagradspolynom
x
f f 2 f 3 f 0 1.80
1.03
0.5
2.83
0.33
1.36
0.11
1.0
4.19
0.44
1.80
1.5
5.99
x f
0 0
ph x h f
1
x
1
Numeriska beräkningar i Naturvetenskap och Teknik
Kvadratisk interpolation Ansats
Q
2 (
x
)
c
0
c
1 (
x
x
0 )
c
2 (
x
x
0 )(
x
x
1 )
1
x
x
0
f
f
0
f
0
c
0
2
x
x
1
f
f
1
f
1
c
0
c
1 (
x
1
x
0 )
c
1
f
1 (
x
1
f
0
x
0 )
f h
3
x
x
2
f
f
2
f
2
c
0
c
1 (
x
2
x
0 )
c
2 (
x
2
x
0 )(
x
2
x
1 )
c
2 (
x
2
f
2
x
0 )(
f x
2 0
x
1 ) (
x f
1 1
f
0
x
0 ) (
x
2 (
x
2
x
0 )(
x
0 )
x
2
x
1 )
Numeriska beräkningar i Naturvetenskap och Teknik
Kvadratisk interpolation 3
c
2 (
x
2
f
2
f
0
x
0 )(
x
2
x
1 )
f
(
x
1 1
f
0
x
0 ) (
x
2 (
x
2
x
0 )(
x
0 )
x
2
x
1 ) (
f
2
f
0 )
h
( 2
h
3
f
1
f
0 ) 2
h
f
2
f
0 2
f
1 2
h
2 2
f
0
f
2 2
f
1 2
h
2
f
0 2
f
2
h
2
Q
2 (
x
)
f
0
p
f
0
f h f
2
f
2
h
(
x
x
0 ) 2
f
2
h
2
p
(
x
x
0
h
) (
x
x
0 )(
x
x
1 )
f
0
p
f
p
(
p
1 ) 2
f
2
h
Numeriska beräkningar i Naturvetenskap och Teknik
Kvadratisk interpolation Newtons ansats
Q m
(
x
)
c
0
c
1 (
x
c
2 (
x
x
0 )
x
0 )(
x
...
x
1 )
c m
(
x
x
0 )(
x
x
1 ) ...
(
x
x m
1 )
f
0
f
1
f
2 ...
f m
Entydighet: Det finns bara ett polynom av gradtal m som går genom m+1 punkter.
Numeriska beräkningar i Naturvetenskap och Teknik
Felet vid interpolation
f
(
x
)
Q
(
x
)
fel
(
x
)
fel
(
x
)
f
(
x
)
Q m
(
x
)
f
(
m
1 ) ( (
m
1 )!
) (
x
x
0 )(
x
x
1 )...(
x
x m
)
Linjär interpolation
fel
(
x
) | |
f
(
x
)
Q
1 (
x
) |
f
'' ( ) | | (
x
2
x
0 )(
x
x
1 ) | | 2
f
2
f
'' ( ) | 2
ph
( 1
p
)
h
h
2 |
p
( 1
p
) (2p 1 0 p 1/2 1/4) 2
f
8
f
2 '' ( ) |
p
( 1
p
)
Numeriska beräkningar i Naturvetenskap och Teknik
Exempel
f
(
x
)
e x
/( 1 9
x
)
Interpolation av polynom av gradtal 4,8,16 i ekvidistanta punkter Anpassning av polynom av gradtal 6 till 9 ekvidistanta punkter
Numeriska beräkningar i Naturvetenskap och Teknik
4 grad 8 grad 16 grad 6 grad i 9 ptr
Numeriska beräkningar i Naturvetenskap och Teknik
Runges fenomen Interpolation i ekvidistanta punkter med ett polynom av högt gradtal tenderar att återge en eftersökt kurva bättre inom de centrala delarna av anpassningsintervallet men ger avsevärda svängningar nära intervallets ändpunkter!
Chebychevpolynom och Chebychevabskissor Om man kan välja de punker i vilka data är kända (kan vara svårt i en given mätserie) så bör mätpunkterna väljas tätare i intervallets ändpunkter. Ett optimalt val ges av Chebychev polynomens nollställen som minimerar resttermen ovan
x i
cos( 2
i m
1 1 2 )
Numeriska beräkningar i Naturvetenskap och Teknik
Splines Ett alternativ är att använda ett polynom av lägre ordning styckvis mellan mätpunkterna. Man kan då ex.vis sätta villkoret att funktions värdena, derivatan och andr derivatan är lika i intervallens ändpunkter för polynom som möts i dessa punkter. Detta ger upphov till s.k. kubiska splines. I ändpunkterna tillkommer kraven ex. vis att kurvan är rak, dvs har konstant derivata utanför intervallet.
f x
Numeriska beräkningar i Naturvetenskap och Teknik
Kubiska splines Funktion
t
: [ 0 , 1 ]
Y i Y i
(
t
) ( 0 )
a y i i
b i t a i
Y i
( 1 )
y i
1
a i c i t
b i
2
d
c i i t
3
d i
Derivatan
Y i Y i
' (
t
) ' ( 0 )
b i D i
2
c i t b i Y i
' ( 1 )
D i
1
b i
3
d i t
2
c i
2 3
d i
Andraderivatan
a b c d i i i i
y i
D i
3 ( 2 (
y i
1
y i
y i
)
y i
1 ) 2
D i D i
D i D i
1 1
Y i
' ' (
t
) 2
c i
6
d i t
Numeriska beräkningar i Naturvetenskap och Teknik
Insättning ger:
Y i
(
t
)
y i
[ 2 (
y i D i t
[ 3 (
y i
1
y i
1 )
D i
y i
) 2
D i D i
1 ]
t
3
D i
1 ]
t
2
Y i
' (
t
)
D i
3 [ 2 ( 2 [ 3 (
y i
1
y i
y i
1 )
y i D i
) 2
D i
D i
1 ]
t
D i
1 ]
t
2
Y i
' ' (
t
) 2 [ 3 (
y i
1
y i
) 2
D i
D i
1 ] 6 [ 2 (
y i
y i
1 )
D i
D i
1 ]
t
Numeriska beräkningar i Naturvetenskap och Teknik
Villkoren
Y i
1
Y i
' 1 ( 1 ) ( 1 )
y i Y
'
i Y i
'' 1 ( 1 )
Y
''
i
( 0 ) ( 0 )
samt
Y Y
'' 0 ( 0 ) 0 ''
n
1 ( 1 ) 0
ger tre okända i tre ekvationer:
y i
y i
1
D i
1 [ 3 (
y i
y i
1 ) 2
D i
1 [ 2 (
y i
1
y i
)
D i
1
D i
]
D i
]
D i
D
i
3 1 [ 2 ( 2 [ 3 (
y i y i
1
y i
)
y i
1 )
D i
1 2
D i
1
D i
]
D i
] 2 [ 3 (
y i
1 2 [ 3 (
y i
y i y i
1 ) ) 2
D i
2
D i
1
D i
1 ]
D i
] 6 [ 2 (
y i
1
y i
)
D i
1
D i
]
Numeriska beräkningar i Naturvetenskap och Teknik
Skrives detta system i matrisform erh:
2 1 1 4 1 1 4 1 ....
1 4 1 2 1 D 0 D 1 D 2 D n 1 D n 3(y 1 3(y 2 3(y 3
y
0 )
y y
0 0 ) ) ....
3 (
y n
3 (
y n
y n
2 )
y n
1 )
Löses enkelt! Testa MATLABs spline funktion själva!