Richardsons extrapolering

Ett exempel på Richardsons extrapolationsmetod i två dimensioner.

I numerisk analys är Richardson-extrapolering en sekvensaccelerationsmetod $A^{\ast }=\lim _{ h\till 0}A(h)$ används för att förbättra konvergenshastigheten för en sekvens av uppskattningar av något värde h . I huvudsak, givet värdet på $A(h)$ för flera värden på $h$ , kan vi uppskatta $A^{\ast }$ genom att extrapolera uppskattningarna till $h=0$ . Den är uppkallad efter Lewis Fry Richardson , som introducerade tekniken i början av 1900-talet, även om idén redan var känd för Christiaan Huygens i hans beräkning av π . Med Birkhoffs och Rotas ord , "kan dess användbarhet för praktiska beräkningar knappast överskattas."

Praktiska tillämpningar av Richardson-extrapolering inkluderar Romberg-integration , som tillämpar Richardson-extrapolation på trapetsregeln , och Bulirsch-Stoer-algoritmen för att lösa vanliga differentialekvationer.

Allmän formel

Notation

Låt $A_{0}(h)$ vara en approximation av $A^{*}$ (exakt värde) som beror på en positiv stegstorlek h med en felformel i formen

$A^{*}=A_{0}(h)+a_{0}h^{k_ {0}}+a_{1}h^{k_{1}}+a_{2}h^{k_{2}}+\cdots$

där $a_{i}$ är okända konstanter och $k_{i}$ är kända konstanter så att $h^{k_{i }}>h^{k_{i+1}}$ . Dessutom representerar $O(h^{k_{i}})$ trunkeringsfelet för $A_{i}(h)$ approximationen så att ${\displaystyle A^{*}=A_{i}(h)+O(h^{k_{i}}).} På liknande sätt,$ i ${ \displaystyle A^{*}=A_{i}(h)+O(h^{k_{i}}),} approximationen$ A $\displaystyle A_{i}(h)}$ sägs vara vara en $O(h^{k_{i}})$ approximation.

Observera att genom att förenkla med Big O-notation är följande formler ekvivalenta:

$${\begin{aligned}A^{*}&=A_{0}(h)+a_{0}h^{k_{0}}+a_{1}h^{ k_{1}}+a_{2}h^{k_{2}}+\cdots \\A^{*}&=A_{0}(h)+a_{0}h^{k_{0}} +O(h^{k_{1}})\\A^{*}&=A_{0}(h)+O(h^{k_{0}})\end{aligned}}$$

Syfte

Richardsons extrapolering är en process som hittar en bättre approximation av $A^{*}$ genom att ändra felformeln från $A^{*}= A_{0}(h)+O(h^{k_{0}})$ till $A^{*}=A_{1}(h)+O(h^{k_{1}}).$ Därför, genom att ersätta $A_{0}(h)$ med $A_{1}(h)$ trunkeringsfelet har minskat från ${\displaystyle O(h^{k_{0}})$ } till ${ \displaystyle O(h^{k_{1}})}$ för samma stegstorlek $h$ . Det allmänna mönstret uppstår där $A_{i}(h)$ har en mer exakt uppskattning än $A_{j}(h)$ när ${ \displaystyle i>j}$ . Genom denna process har vi uppnått en bättre approximation av $A^{*}$ genom att subtrahera den största termen i felet som var $O(h^{k_{0}} )$ . Denna process kan upprepas för att ta bort fler feltermer för att få ännu bättre uppskattningar.

Bearbeta

Genom att använda stegstorlekarna $h$ och $h/t$ för någon konstant ${\displaystyle t} ,$ är de två formlerna för ${\displaystyle A^{*}} :$

$${\begin{aligned}A^{*}&=A_{0}(h)+ a_{0}h^{k_{0}}+a_{1}h^{k_{1}}+a_{2}h^{k_{2}}+O(h^{k_{3}}) &(1)\\\\A^{*}&=A_{0}\!\left({\frac {h}{t}}\right)+a_{0}\left({\frac {h }{t}}\right)^{k_{0}}+a_{1}\left({\frac {h}{t}}\right)^{k_{1}}+a_{2}\left ({\frac {h}{t}}\right)^{k_{2}}+O(h^{k_{3}})&(2)\end{aligned}}$$

För att förbättra vår approximation från $O(h^{k_{0}})$ till $O(h^{k_{1}})$ genom att ta bort första feltermen multiplicerar vi den andra ekvationen (2) med $t^{k_{0}}$ och subtraherar den första ekvationen (1) för att ge oss

$$(t^{k_{0}}-1)A^{*}={\bigg [}t^{k_{0}}A_{0}\left({\frac {h}{t} }\right)-A_{0}(h){\bigg ]}+{\bigg (}t^{k_{0}}a_{1}{\bigg (}{\frac {h}{t}} {\bigg )}^{k_{1}}-a_{1}h^{k_{1}}{\bigg )}+{\bigg (}t^{k_{0}}a_{2}{\ bigg (}{\frac {h}{t}}{\bigg )}^{k_{2}}-a_{2}h^{k_{2}}{\bigg )}+O(h^{k_ {3}}).$$

Denna multiplikation och subtraktion utfördes för att göra ${\textstyle {\big [}t^{k_{0}}A_{0}\left({\frac {h }{t}}\right)-A_{0}(h){\big ]}}$ en $O(h^{k_{1}})$ approximation av $(t^{k_{0}}-1)A^{*}$ . Vi kan lösa vår nuvarande formel för $A^{*}$ att ge

$$A^{*}={\frac {{\bigg [}t^{k_{0}} A_{0}\left({\frac {h}{t}}\right)-A_{0}(h){\bigg ]}}{t^{k_{0}}-1}}+{\ frac {{\bigg (}t^{k_{0}}a_{1}{\bigg (}{\frac {h}{t}}{\bigg )}^{k_{1}}-a_{1 }h^{k_{1}}{\bigg )}}{t^{k_{0}}-1}}+{\frac {{\bigg (}t^{k_{0}}a_{2} {\bigg (}{\frac {h}{t}}{\bigg )}^{k_{2}}-a_{2}h^{k_{2}}{\bigg )}}{t^{ k_{0}}-1}}+O(h^{k_{3}})$$

som kan skrivas som $A^{*}=A_{1}(h)+O(k^{k_{1}})$ av inställning

$A_{1}(h)={\frac {t^{k_{0}}A_{0}\left({\frac {h}{t}}\right)-A_{0}(h )}{t^{k_{0}}-1}}.$

Återkommande förhållande

En allmän återkommande relation kan definieras för approximationerna med

$A_{i+1}(h)={\frac {t^{k_{ i}}A_{i}\left({\frac {h}{t}}\right)-A_{i}(h)}{t^{k_{i}}-1}}$

där $k_{i+1}$ uppfyller

$A^{*}=A_{i+1}(h)+O(h^{k_{i+1}} )$ .

Egenskaper

Richardson-extrapoleringen kan betraktas som en linjär sekvenstransformation .

Dessutom kan den allmänna formeln användas för att uppskatta $k_{0}$ (stegstorleksbeteende i ledande ordning för Truncation error ) när varken dess värde eller $A^{*}$ är känt a priori . En sådan teknik kan vara användbar för att kvantifiera en okänd konvergenshastighet . Givet approximationer av $A^{*}$ från tre distinkta stegstorlekar $h$ , $h/t$ , och $h/s$ , exakt förhållande

$$A^{*}={\frac {t^{k_{0}}A_{i}\left({\frac {h}{t}}\right)-A_{i} (h)}{t^{k_{0}}-1}}+O(h^{k_{1}})={\frac {s^{k_{0}}A_{i}\left({ \frac {h}{s}}\right)-A_{i}(h)}{s^{k_{0}}-1}}+O(h^{k_{1}})$$

ger ett ungefärligt samband (observera att notationen här kan orsaka lite förvirring, de två O som visas i ekvationen ovan indikerar bara beteendet för stegstorleken i ledande ordning, men deras explicita former är olika och därför är det att ta bort de två O-termerna ungefär giltigt)

$$A_{i}\left({\frac {h}{t}}\right)+{\frac {A_{i}\left({\frac {h}{t}}\right)-A_{i} (h)}{t^{k_{0}}-1}}\approx A_{i}\left({\frac {h}{s}}\right)+{\frac {A_{i}\left ({\frac {h}{s}}\right)-A_{i}(h)}{s^{k_{0}}-1}}$$

som kan lösas numeriskt för att uppskatta $k_{0}$ för några godtyckliga val av $h$ , $s$ och $t$ .

Exempel på Richardson-extrapolering

Antag att vi vill approximera $A^{*}$ , och vi har en metod $A(h)$ som beror på en liten parameter $h$ i en sådan sätt att

$$A(h)=A^{\ast }+Ch^{n}+O(h^{n+1}).$$

Låt oss definiera en ny funktion

$$R(h,t):={\frac {t^{n}A(h/t )-A(h)}{t^{n}-1}}$$

där $h$ och ${\frac {h}{t}}$ är två distinkta stegstorlekar.

Sedan

$$R(h,t)={\frac {t^{n}(A^{*}+C\left({\frac {h}{t}}\right)^{n}+O( h^{n+1}))-(A^{*}+Ch^{n}+O(h^{n+1}))}{t^{n}-1}}=A^{* }+O(h^{n+1}).$$

$R(h,t)$ kallas Richardson- extrapoleringen av A ( h ), och har en högre ordningens feluppskattning $O(h^{ n+1})$ jämfört med $A(h)$ .

Mycket ofta är det mycket lättare att få en given precision genom att använda R ( h ) snarare än A ( h′ ) med mycket mindre h′ . Där A ( h′ ) kan orsaka problem på grund av begränsad precision ( avrundningsfel ) och/eller på grund av det ökande antalet beräkningar som behövs (se exempel nedan).

Exempel pseudokod kod för Richardson extrapolation

Följande pseudokod i MATLAB-stil demonstrerar Richardsons extrapolering för att hjälpa till att lösa ODE $y'(t)=-y^{2}$ , $y( 0)=1$ med trapetsmetoden . I det här exemplet halverar vi stegstorleken $h$ varje iteration och så i diskussionen ovan skulle vi ha att $t=2$ . Felet i den trapetsformade metoden kan uttryckas i termer av udda potenser så att felet över flera steg kan uttryckas i jämna potenser; detta leder till att vi höjer $t$ till andra potensen och tar potenserna $4=2^{2}=t^{2}$ i pseudokoden. Vi vill hitta värdet på $y(5)$ , som har den exakta lösningen $​​{\displaystyle {\frac {1}{5+ 1}}={\frac {1}{6}}=0,1666...}$ eftersom den exakta lösningen av ODE är $y(t)={\frac {1 }{1+t}}$ . Denna pseudokod antar att en funktion som kallas Trapezoidal(f, tStart, tEnd, h, y0) existerar som försöker beräkna y(tEnd) genom att utföra trapetsmetoden på funktionen f , med startpunkt y0 och tStart och stegstorlek h .

Observera att att börja med en för liten inledande stegstorlek potentiellt kan introducera fel i den slutliga lösningen. Även om det finns metoder utformade för att hjälpa till att välja den bästa initiala stegstorleken, är ett alternativ att börja med en stor stegstorlek och sedan låta Richardson-extrapoleringen minska stegstorleken varje iteration tills felet når den önskade toleransen.

  0          
              
              
                    
                
    

                  
        
   

  




   


       



        
                   
    
    
             
    
              
        
        
        
     
                       
    
    
    
    
    
    
                 
                                       
          
                                                          
    


      
     
      
 tStart  =  % Starttid  tEnd  =  5  % Sluttid  f  =  -  y  ^  2  % Derivatan av y, så y' = f(t, y(t)) = -y^2  % Lösningen till denna ODE är y = 1/(1 + t)  y0  =  1  % Startpositionen (dvs. y0 = y(tStart) = y(0) = 1)  tolerans  =  10  ^  -  11  % 10-siffrig noggrannhet önskas  max . Rader   =  20  % Tillåt inte iterationen för att fortsätta på obestämd tid  initialH  =  tStart  -  TEnd  % Välj ett inledande stegstorlek  haveWeFoundSolution  =  falskt  % Kunde vi hitta lösningen inom den önskade toleransen? inte än.   h  =  initialH  % Skapa en 2D-matris med storleken maxRows by maxRows för att hålla Richardson-extrapoleringen  % Observera att detta kommer att vara en lägre triangulär matris och att högst två rader faktiskt behövs  % när som helst i beräkningen.  A  =  nollMatrix  (  maxRows  ,  maxRows  )  % Beräkna det övre vänstra elementet i matrisen. Den första raden i denna (nedre triangulära) matris har nu fyllts.   A  (  1  ,  1  )  =  trapetsformad  (  f  ,  tStart  ,  tEnd  ,  h  ,  y0  )  % Varje rad i matrisen kräver ett anrop till trapetsformad %  Denna loop börjar med att fylla den andra raden i matrisen, eftersom den första raden beräknades ovan  för  i  =  1  :  maxRows  -  1  % Börja vid i = 1, iterera högst maxRows - 1 gånger  h  =  h  /  2  % Halvera det föregående värdet på h eftersom detta är början på en ny rad.  % Börja fylla rad i+1 från vänster genom att anropa den trapetsformade funktionen med denna nya mindre stegstorlek  A  (  i  +  1  ,  1  )  =  Trapetsformad  (  f  ,  tStart  ,  tEnd  ,  h  ,  y0  )  för  j  =  1  :  i  % Go över denna nuvarande (i+1)-te rad tills diagonalen nås  % För att beräkna A(i + 1, j + 1), vilket är nästa Richardson-extrapolering,  använder % det senast beräknade värdet (dvs. A(i + 1, j)) och värdet från  %-raden ovanför den (dvs. A(i, j)).  A  (  i  +  1  ,  j  +  1  )  =  ((  4  ^  j  )  .*  A  (  i  +  1  ,  j  )  -  A  (  i  ,  j  ))  /  (  4  ^  j  -  1  );  slut  % Efter att ha lämnat ovanstående inre slinga har diagonalelementet i rad i + 1 beräknats  % Detta diagonalelement är det senaste Richardson-extrapoleringen som har beräknats.  % Skillnaden mellan denna extrapolering och den sista extrapoleringen av rad i är en bra  % indikation på felet.  if  (  absolutVärde  (  A  (  i  +  1  ,  i  +  1  )  -  A  (  i  ,  i  ))  <  tolerans  )  % Om resultatet är inom tolerans  skriv ut  (  "y = "  ,  A  (  i  +  1  ,  i  +  1  ))  % Visa resultatet av Richardson-extrapoleringen  haveWeFoundSolution  =  true  break  %  Klart, så lämna slingslutet  if  (  haveWeFoundSolution  ==  false  )  %  Om vi ​​inte kunde hitta en lösning inom den önskade toleransutskriften  (  "  Warning: Not able att hitta en lösning inom den önskade toleransen av " ,  tolerans  )  ;  print  (  "Den sista beräknade extrapoleringen var "  ,  A  (  maxRows  ,  maxRows  ))  slut 

Se även

• Aitkens delta-kvadratprocess
• Takebe Kenko
• Richardson iteration

• Extrapolationsmetoder. Teori och praktik av C. Brezinski och M. Redivo Zaglia , North-Holland, 1991.
• Ivan Dimov, Zahari Zlatev, Istvan Farago, Agnes Havasi: Richardson Extrapolation: Practical Aspects and Applications ,   Walter de Gruyter GmbH & Co KG, ISBN 9783110533002 (2017).

externa länkar

• Grundläggande metoder för numerisk extrapolering med applikationer, mit.edu
• Richardson-Extrapolation
• Richardsons extrapolering på en webbplats för Robert Israel (University of British Columbia)
• Matlab-kod för generisk Richardson-extrapolering.
• Julia-kod för generisk Richardson-extrapolering.