Demo 5: Taylor-approksimeringer

Demo af Christian Mikkelstrup, Hans Henrik Hermansen, Jakob Lemvig, Karl Johan Måstrup Kristensen og Magnus Troen. Revideret marts 2026 af shsp.

from sympy import *
from sympy.abc import x,y,z,u,v,w,t
from dtumathtools import *
init_printing()

Taylor-polynomier af én variabel

I Sympy udføres Taylor-udvikling med kommandoen series(funktion, variabel, udviklingspunkt, grad+1) eller tilsvarende funktion.series(variabel, udviklingspunkt, grad+1). En approksimering af \(\cos(x)\) med et 6.-grads Taylor-polynomium fra udviklingspunktet \(x_0 = 0\) fås med koden:

series(cos(x), x, 0, 7)

Bemærk, hvordan et approksimerende polynomium af grad 6 opnås ved at taste 7 som det sidste argument i series-kommandoen. Generelt set skal man taste den ønskede grad plus 1 i kommandoen, hvilket kan føles spøjst men er vigtigt at huske. Den tekniske årsag er, at outputtet, som du ser herover, ikke kun indeholder det approksimerende polynomium men også et restled, der i den viste form altid vil være af én højere orden end polynomiets grad. Det er restleddets orden, series-kommandoen kræver som argument.

“Store O”-notation

Restleddet vist i outputtet ovenfor er formentlig ikke noget, du har set før. Lærebogen præsenterer restled enten som epsilon-funktioner eller som (Laplace)-restfunktioner, der er en del af kursets pensum, men her ser vi såkaldt “Store O”-notation, der for Taylor-polynomier af grad \(K\) har formen \(O(x^{K+1})\). Opfat dette som blot endnu en fejlindikator men med et andet formål end de andre.

Hvor fx restfunktionen kan benyttes til fejlvurdering, fortæller “Store O” noget om vækstraten af fejlen. Dvs., “Store O” giver en øvre grænse for, hvor hurtigt fejlen vil kunne vokse, hvis du bevæger dig længere væk fra udviklingspunktet. Outputtet ovenfor fortæller, at fejlen mellem det approksimerende polynomiums værdi og den originale, ægte funktionsværdi med \(\cos(x)\) forbliver under \(x^7\). Alternativt: Som du nærmer dig udviklingspunktet, formindskes fejlen hurtigere end \(x^7\) formindskes.

I dette kursus arbejder vi ikke med vækstrater for den indførte fejl, så hvis kun selve Taylor-polynomiet uden restledsinformation ønskes, kan \(O(x^{K+1})\)-leddet fjernes med .removeO():

series(cos(x), x, 0, 7).removeO()

Eksempel med sin(x)

Betragt funktionen \(f(x)=\sin(x)\) approksimeret fra udviklingspunktet \(x_0 = 0\).

Lad os plotte \(f(x)\) sammen med nogle udvalgte approksimationer. Vi benytter argumentet show = False i plotkommandoen, så vi kan tilføje flere plots til plotvinduet med .extend(), før det samlede resultat vises:

pl = plot(sin(x),xlim = (-3,3), ylim = (-3,3), show=False, legend = True)
for K in [0,1,2,3,6]:
    nyraekke = series(sin(x),x,0,K+1).removeO()
    display(Eq(Function(f'P_{K}')(x), nyraekke))
    nytplot = plot(nyraekke,label = f"n = {K}", show=False)
    pl.extend(nytplot)
pl.show()

Det kombinerede plot viser tydeligt, at en forøgelse af graden \(K\) giver en bedre og bedre approksimation. Bemærk, at nogle approksimationer overlapper:

\[\begin{align*} P_1(x) &= f(x_0) + f'(x_0)(x-x_0)\\ &= \sin(0) + \cos(0) x \\ &= x \\ P_2(x) &= f(x_0) + f'(x_0)(x-x_0) + \frac{1}{2} f''(x_0)(x-x_0)^2 \\ &= \sin(0) + \cos(0) x - \sin(0) x^2 \\ &= x \end{align*}\]

for \(x\in\mathbb{R}\).

Fejlvurdering

Vi ønsker at finde funktionsværdien \(\ln\left(\frac{5}{4}\right)\), men den kan være besværlig at beregne nøjagtigt (især hvis vi skal finde mange af sådanne værdier). Lad os gøre det nemmere ved i stedet at benytte os af et approksimerende tredjegradspolynomium \(P_3(x)\) udviklet fra punktet \(x_0=1\):

x = symbols("x")
P3 = series(ln(x),x,1,4).removeO()
P3

Vi kan nu bruge denne til at finde en approksimerende værdi, altså en “cirkaværdi”, af \(\ln(x)\) i \(x=\frac{5}{4}\):

vaerdi = P3.subs(x,Rational(5/4))
vaerdi, vaerdi.evalf()

Altså, \(\ln\left(\frac{5}{4}\right)\approx 0.224\). Men hvor nøjagtig er denne værdi? Hvor stor en fejl har vi indført ved at benytte approksimationen i stedet for den ægte funktion? Det vil vi nu undersøge.

Maksimering af restfunktionen

Lærebogen fortæller, at et \(\xi \in \,\,]\,x_0\,,\,x\,[\) eksisterer, således at fejlen kan skrives som en restfunktion \(R_n(x)\):

\[\begin{equation*} R_n\left(x\right) = \frac{f^{(n+1)}(\xi)}{(n+1)!}\cdot\left(x - x_0\right)^{n+1}. \end{equation*}\]

Denne restfunktion repræsenterer den indførte fejl, som vi begår ved at anvende approksimationen i stedet for den sande funktion. Pakket ind i absolutværdistreger, \(|R_n(x)|\), repræsenterer den den numeriske fejl (uden fortegn). Den største værdi, som (den absolutte værdi af) restfunktionen kan antage, er dermed en øvre grænse for fejlen. Den faktiske, nøjagtige fejl, der begås, er muligvis mindre, men da vi ikke er i stand til at finde den nøjagtige fejlværdi, er det meget brugbart at have information om den øvre grænse for fejlen.

Lad os indsætte, hvad vi ved: I \(x=\frac{5}{4}\) med \(x_0=1\) bliver intervallet for \(\xi\) til \(\xi \in \,\,]\,1\,,\,\frac{5}{4}\,[\), og fejlfunktionen bliver:

\[\begin{equation*} R_3\left(\frac{5}{4}\right) = \frac{f^{(4)}(\xi)}{4!}\cdot\left(\frac{5}{4} - 1\right)^4. \end{equation*}\]

Lad os udregne \(f^{(4)}(\xi)\):

xi = symbols("xi")
diff(ln(x),x,4).subs(x,xi)

Vi har dermed:

\[\begin{equation*} R_3\left(\frac{5}{4}\right) = \frac{-6}{4!\xi^4}\cdot\left(\frac{5}{4} - 1\right)^4. \end{equation*}\]

Den eneste ubekendte er \(\xi\), og den er tilfældigvis i nævneren her. Derfor maksimeres \(\left|R_3\left(\frac{5}{4}\right)\right|\) ved en minimering af \(\xi\). I sit interval vælger vi derfor den mindst mulige værdi, så \(\xi=1\), som indsættes:

R3 = abs(diff(ln(x),x,4).subs(x,1) * (5/4 - 1) ** 4 /(factorial(4)))
Eq(Function('R_3')(S('5/4')), R3)

Vi har nu en øvre grænse for den numeriske fejl. Den sande (numeriske) fejl vil være mindre end dette. Hvis nødvendigt kan vi opstille et fejlinterval med endepunkterne:

Interval(vaerdi - R3, vaerdi + R3)

Altså, konklusioner er: Vi kan garantere, at den nøjagtige værdi af \(\ln(\frac{5}{4})\) ligger et sted i intervallet \(]0.2229;0.2250[\) (her et åbent interval grundet afrundning). Det kunne nu være interessant at tjekke, om Python er enig - vi beder Python om at udregne den sande værdi af \(\ln(\frac{5}{4})\) for at bekræfte, at den ganske rigtigt er indenfor den fundne fejlmargin fra den approksimerende værdi:

ln(5/4), Interval(vaerdi - R3, vaerdi + R3).contains(ln(5/4))

Bemærk, at selv Pythons værdi for \(\ln(\frac{5}{4})\) er en approksimation, blot en langt bedre en end vores. Vi kan dog sagtens forbedre vores approksimation ved simpelthen at udføre højereordens Taylor-approksimationer.

Grænseværdier med Taylors grænseformel

Lad os nu anvende Taylors grænseformel til at bestemme grænseværdier af forskellige udtryk. Dette er ofte brugbart til brøker, hvori tæller eller nævner (eller begge) indeholder udtryk, der ikke er nemme at arbejde med. For eksempel er det ikke let at se, hvilken værdi udtrykket \(\frac{\sin(x)}{x}\) går imod for \(x\rightarrow 0\), da både tæller og nævner indeholder \(x\).

Vi tilgår Taylors grænseformel i Sympy med series-kommandoen, der denne gang ikke skal tilføjes .removeO(), da vi nu har brug for et restudtryk. Bemærk, at kursets lærebog præsenterer Taylors grænseformel med en epsilonfunktion i restleddet, hvorimod Sympy bruger “Store O”-notationen som omtalt tidligere. Vi konverterer Sympys output til epsilonformatet manuelt ved at erstatte “Store O”-notation som \(O(x^{K+1})\) med \(\varepsilon(x) \, x^{K}\), hvor \(\varepsilon(x) \to 0\) for \(x \to 0\).

Eksempel 1

Vi vil undersøge grænseværdien af udtrykket \(\frac{\sin(x)}{x}\) for \(x\rightarrow 0\). Lad os Taylor-udvikle tælleren til anden orden:

series(sin(x),x,0,n=3)

Husk på at en sådan Taylor-udvikling er identisk med den oprindelige funktion, ikke blot en approksimation af den (netop fordi restleddet er inkluderet). Vi kan derfor erstatte funktionen direkte med sin Taylor-udvikling, hvilket pludselig lader os reducere brøken og finde dens grænseværdi:

\[\begin{gather*} \frac{\sin(x)}{x} = \frac{x + \epsilon (x) \cdot x^2}{x} = 1 + \varepsilon (x)x \rightarrow 1, \text{ når } x \rightarrow 0. \end{gather*}\]

Dette virker, fordi vi startede med polynomieled i nævneren, og via Taylor-udviklingen konverterede vi også tælleren til polynomieled (plus en epsilonfunktion). Leddene har da mulighed for at gå ud med hinanden, så brøken forhåbentlig kunne reduceres ned til en form, hvor vi ikke længere havde et \(x\) i både tæller og nævner.

Eksempel 2

Hvordan vidste vi, hvilken orden der skulle Taylor-udvikles til i forrige eksempel? Det føltes måske som gætteværk at finde lige netop det rette antal af lige netop de rette led.

Betragt dette udtryk:

\[\begin{equation*} \frac{\mathrm e^x - \mathrm e^{-x} - 2x}{x-\sin(x)}. \end{equation*}\]

Vi ønsker at finde dets grænseværdi for \(x \rightarrow 0\). Denne gang indeholder både tæller og nævner mere end bare polynomieled, så lad os Taylor-udvikle begge hver for sig. Vi vil da få at se, hvad der sker, som vi øger ordenerne og derved tilfører flere og flere led. Vi lægger ud med at udvikle begge til 1. orden:

T = E ** x - E ** (-x) - 2*x
N = x - sin(x)
series(T,x,0,2), series(N,x,0,2)

Dette er for upræcist, da Taylor-udviklingerne blot bliver restleddene. Indsættelse af disse i brøken giver os ikke mulighed for at reducere yderligere (epsilonfunktionerne for hver Taylor-udvikling kan være forskellige, så de går ikke ud med hinanden). Lad os prøve med 2. orden:

series(T,x,0,3), series(N,x,0,3)

Stadig for upræcist. Hvad så med til 3. orden:

series(T,x,0,4), series(N,x,0,4)

Her har vi noget, vi kan forsøge med. Indsættelse i brøken giver:

\[\begin{gather*} \frac{\frac{x^3}{3}+\epsilon_1(x)x^3}{\frac{x^3}{6}+\epsilon_2(x)x^3} = \frac{\frac{1}{3} + \epsilon_1(x)}{\frac{1}{6}+\epsilon_2(x)} \rightarrow \frac{\frac{1}{3}}{\frac{1}{6}} = 2, \text{ når } x \rightarrow 0. \end{gather*}\]

Vi har sat et indeks på hver epsilonfunktion, da de kan være forskellige. Men det har ingen betydning, da de begge går imod nul ved grænsen. Lad os kontrollere med Pythons limit()-kommando:

limit(T/N,x,0)

Taylor-polynomier af to variable

Betragt følgende funktion af to variabler:

\[ f:\mathbb{R}^2 \to \mathbb{R},\quad f(x,y) = \sin(x^2 + y^2). \]

Lad os først og fremmest plotte den:

x,y = symbols("x y", real = True)
f = sin(x ** 2 + y ** 2)
dtuplot.plot3d(f,(x,-1.5,1.5),(y,-1.5,1.5),rendering_kw={"color" : "blue"})

<spb.backends.matplotlib.matplotlib.MatplotlibBackend at 0x7f183aab9e10>

Approksimerende polynomier til forskellige ordener

Lad os finde funktionens approksimerende førstegradspolynomium \(P_1(x,y)\) fra udviklingspunkt \((0,0)\):

P1 = dtutools.taylor(f,[x,0,y,0],degree=2)
P1

\(P_1(x,y)\) viser sig blot at være nul, altså et horisontalt plan i højden nul (identisk med \((x,y)\)-planet):

p = dtuplot.plot3d(P1,(x,-1.5,1.5),(y,-1.5,1.5),show=false,rendering_kw={"alpha" : 0.5},camera={"azim":-81,"elev":15})
p.extend(dtuplot.plot3d(f,(x,-1.5,1.5),(y,-1.5,1.5),show=False))
p.show()

Med et andet udviklingspunkt så som \((1/10,0)\) får vi:

p = dtuplot.plot3d(dtutools.taylor(f,[x,0.1,y,0],degree=2),(x,-1.5,1.5),(y,-1.5,1.5),show=false,rendering_kw={"alpha" : 0.5},camera={"azim":-81,"elev":15})
p.extend(dtuplot.plot3d(f,(x,-1.5,1.5),(y,-1.5,1.5),show=False))
p.show()

hvilket tydeligt viser, at udviklingspunktets placering har stor betydning.

Lad os vende tilbage til udviklingspunktet \((0,0)\) og undersøge det approksimerende andengradspolynomium:

P2 = dtutools.taylor(f,[x,0,y,0],3)
P2

dtuplot.plot3d(f,P2,(x,-1.5,1.5),(y,-1.5,1.5))

<spb.backends.matplotlib.matplotlib.MatplotlibBackend at 0x7f1833e76a10>

Denne gang er det approksimerende polynomium en (elliptisk) paraboloide.

Lad os slutteligt gå en smule højere op i orden og tage et kig på det approksimerende sjettegradspolynomium:

P6 = dtutools.taylor(f,[x,0,y,0],7)
p = dtuplot.plot3d(f,P6,(x,-1.5,1.5),(y,-1.5,1.5),show=False)
p.legend=True
p.show()

Som forventet stemmer dette langt bedre med funktionen.

Fejlundersøgelse

For at se, hvor godt de stemmer, vil vi undersøge fejlene, vi begår, når vi bruger et approksimerende polynomium af grad 1, 2 og 6 i stedet for den originale funktion til at evaluere funktionsværdier i forskellige punkter. Vi lægger ud med punktet \((0.2,0.2)\):

f_p1 = f.subs([(x, 1/5), (y, 1/5)])
P1_p1 = P1.subs([(x, 1/5), (y, 1/5)])
P2_p1 = P2.subs([(x, 1/5), (y, 1/5)])
P6_p1 = P6.subs([(x, 1/5), (y, 1/5)])

fejl_liste = (f_p1 - P1_p1, f_p1 - P2_p1, f_p1 - P6_p1)
klargjorte_ligninger = [ Eq(Function(f'P_{i}')(S('1/5')), udtryk) for i, udtryk in zip((1,2,6), fejl_liste) ]

display(*klargjorte_ligninger)

Som forventet er fejlen langt mindre for højereordens approksimationer. Lad os prøve i \((0.5,0.5)\):

f_p2 = f.subs([(x,1/2),(y,1/2)])
P1_p2 = P1.subs([(x,1/2),(y,1/2)])
P2_p2 = P2.subs([(x,1/2),(y,1/2)])
P6_p2 = P6.subs([(x,1/2),(y,1/2)])

fejl_liste = (f_p2 - P1_p2, f_p2 - P2_p2, f_p2 - P6_p2)
klargjorte_ligninger = [ Eq(Function(f'P_{i}')(S('1/5')), udtryk) for i, udtryk in zip((1,2,6), fejl_liste) ]

display(*klargjorte_ligninger)

Jo længere væk fra udviklingspunktet \((0,0)\) vi kommer, desto større bliver fejlen, præcis som forventet.

Bemærkning: Det bør nævnes, at vores sammenligninger her er baseret på en antagelse om, at Sympys egne approksimeringer er langt bedre end vores. Det er højst sandsynligt en ganske fornuftig antagelse i tilfælde som disse, men det er vigtigt at huske på, at Sympy (for så vidt computer software generelt) også udfører approksimationer.

Taylor-polynomier af tre variables

Lad os i dette afsluttende afsnit undlade at bruge series-kommandoen og i stedet forsøge at konstruere en andenordens Taylor-approksimation (semi)manualt. Til det formål giver lærebogen os følgende formel for andengrads Taylor-polynomier af funktioner i flere variable:

\[\begin{equation*} P_2(\boldsymbol{x}) = f(\boldsymbol{x}_0) + \left<(\boldsymbol{x} - \boldsymbol{x}_0), \nabla f(\boldsymbol{x}_0)\right> + \frac{1}{2}\left<(\boldsymbol{x} - \boldsymbol{x}_0), \boldsymbol{H}_f(\boldsymbol{x}_0)(\boldsymbol{x}-\boldsymbol{x}_0)\right>. \end{equation*}\]

Lad os som et eksempel betragte funktionen:

\[\begin{equation*} f: \mathbb{R}^3\to \mathbb{R},\quad f(x_1,x_2,x_3) = \sin(x_1^2 - x_2)\mathrm e^{x_3}. \end{equation*}\]

x1,x2,x3 = symbols('x1,x2,x3', real = True)
f = sin(x1**2 - x2)*exp(x3)
f

Vi ønsker at bestemme andengrads Taylor-polynomiet med udviklingspunkt \(\boldsymbol{x}_0 = (1,1,0)\). Først defineres \(\boldsymbol{x}_0\) og \(\boldsymbol{x}\):

x = Matrix([x1,x2,x3])
x0 = Matrix([1,1,0])

Dernæst udregnes \(\nabla f(\boldsymbol{x}_0)\) og \(\boldsymbol{H}_f(\boldsymbol{x}_0)\):

nabla_f = dtutools.gradient(f,(x1,x2,x3)).subs([(x1,x0[0]),(x2,x0[1]),(x3,x0[2])])
nabla_f

\[\begin{split}\displaystyle \left[\begin{matrix}2\\-1\\0\end{matrix}\right]\end{split}\]

Hf = dtutools.hessian(f,(x1,x2,x3)).subs([(x1,x0[0]),(x2,x0[1]),(x3,x0[2])])
Hf

\[\begin{split}\displaystyle \left[\begin{matrix}2 & 0 & 2\\0 & 0 & -1\\2 & -1 & 0\end{matrix}\right]\end{split}\]

Vi har nu alle nødvendige komponenter til konstruering af polynomiet \(P_2\):

P2 = f.subs([(x1,x0[0]),(x2,x0[1]),(x3,x0[2])]) + nabla_f.dot(x - x0) + S('1/2')* (x - x0).dot(Hf*(x - x0))
P2.simplify()

Lad os igen tage et hurtigt blik på fejlene, som vi begår ved at bruge denne approksimation i stedet for den oprindelige funktion i nogle udvalgte punkter:

v1 = Matrix([1,1,0])
v2 = Matrix([1,0,1])
v3 = Matrix([0,1,1])
v4 = Matrix([1,2,3])
vs = [v1,v2,v3,v4]

for v in vs:
    print((f.subs({x1:v[0],x2:v[1],x3:v[2]}) - P2.subs({x1:v[0],x2:v[1],x3:v[2]})).evalf())

0
0.287355287178842
0.712644712821158
-12.9013965351501

Igen ser vi, at fejlen stiger, som vi bevæger os længere væk fra udviklingspunktet, helt som forventet.