140 89 17MB
Swedish Pages [380] Year 1995
Föreläsningar i Lineära system
Sven Spanne
Matematiska Institutionen
LUIJds Tekniska Högskola
Tryck och distribution: KFS AB Sölvegatan 22 223 62 Lund Tel. 046-184191 Fax 046-140747
Mångfaldigandet av innehållet i denna bok, helt eller delvis, är enligt lagen om upphovsrätt av den 30 december 1960 förbjudet utan medgivande av copyrightinnehavaren.
© Sven Spanne 1982, 1993, 1995 Tredje upplagan 1995 Fjärde tryckningen 1997
Företal
Ur företalet till andra upplagan Innehållet i kursen Lineära system har i stort sett varit oförändrat sedan den gavs för första gången hösten 1981. Under denna tid har utvecklingen av persondat~ rerna inneburit en revolution i möjligheterna för ingenjören att använda matematik. Med numeriska och symboliska program som Matlab och Maple kan man nu på sitt skrivbord enkelt utföra beräkningar som för tjugo år var nästan otänkbara. Den teoretiska bakgrunden är dock fortfarande densamma, och blir allt viktigare, ju mer räkning man överlämnar åt datorerna. Detta gör att kursen inte förlorat något i aktualitet. Möjligheterna att låta datorer sköta beräkningarna gör inte heller lösandet av enkla övningsexempel för hand överflödigt. I denna liksom i många andra matematikkurser på högre nivå är övningarna snarast till för att man skall förstå teorin och begreppen och inte tvärt om. Läsaren uppmanas dock provat ex Matlab och Maple (som finns tillgängliga på skolans datorer för dem som läser kursen vid LTH) och se hur kraftfulla verktyg sådana program är. Kursen läses efter sammanlagt en termins (20 poäng) matematikstudier vid högskolan. De områden som berörs hör i många fall till den högre matematiken, och det är klart att man inte kan vänta sig en djup och fullständig behandling av problemen. Dessutom har läsarna inte den bakgrund i tekniska ämnen, som skulle behövas för att få en full förståelse för problemställningarna. Jag har därför använt mycket utrymme för att diskutera relationerna mellan matematiken och tillämpningarna. Det skall dock påpekas att kursen är ingen lärobok i kretsteori eller systemteknik, utan tillämpningsexemplen skall ses som smakprov på vad matematik kan användas till. Det är också så att alla de metoder som finns beskrivna används lika mycket it ex klassisk och modern fysik som i reglerteknik och signalteori. Några av mina kollegor har bistått med korrekturläsning på ett tidigt stadium. och det är deras förtjänst att ett stort antal tryckfel och oklarheter kunnat avlägnas lll
iv
i tid. Till dem framför jag ett varmt tack. och hoppas att läsarna har glädje av deras insatser.
Lund, i hast
Sven Spanne
Förord till tredje upplagan Den gamla kursen i Lineära system vid LTH har nu gått till historien, men dess innehåll (och anda) lever kvar till stor del i de nya kurserna i Lineär analys. Därför har denna bok liksom fågel Fenix fått ett nytt liv. Ännu har inte multimedia och hypertext slagit helt igenom. och den tredje upplagan skiljer sig från den andra endast genom måttliga förbättringar. tillägg och rättelser. Jag upprepar mitt tack till alla som bidragit genom förslag till korrigeringar. Vid denna revision kan speciellt nämnas Erik Andersen, Tomas Carnstam och Arne Persson.
Lund i höstrusket 1995
Sven Spanne
Innehåll Företal
ili
Innehåll
V
1 Dynamiska system 1. 1 Matematisk beskrivning av system 1.2 Statiska och dynamiska system .. 1.3 Tillståndsvariabler . . . . . . . . . 1.4 Konstruktion av tillståndsmodeller .
1 1 2 3
2 Lineära tillståndsmodeller 2.1 Lineära och olineära funktioner . 2.2 Källor till lineära modeller 2.3 Linearisering . . . . . . . . 2.4 Linearisering . . . . . . . . 2.5 Lineära tillståndsmodeller . 2.6 Derivator av högre ordning 2.7 Kirchhoffs lagar 2. 8 Referenser . . . . . . . . .
9 9
3 Differentialekvationer och spektralteori 3.1 Separation av variabler och superposition 3.1.1 Separation av variabler . . . . 3.2 Egenvärden och egenvektorer . . . . 3.3 Exempel på egenvärdesbestämning . 3.4 Byte av tillståndsvariabler . . . 3.5 Diagonalisering . . . . . . . . . . 3.6 Lösning av tillståndsekvationer . 3. 7 Icke diagonaliserbara matriser . 3.8 Komplexa egenvärden och reella lösningar . Några egenskaper hos egenvärden 3.9 3.10 Fouriers metod . 3.11 Referenser . . . . . . . . . . V
4
12 12 16 18
24 25 27
29 29
32 34 36 39 41 44 49
52 55 57 59
vi
INNEHÅLL
4 Stabilitet och stationära lösningar 4.1 Stabilitet av system 4.2 Fria svängningar . . . 4.3 Tvungna svängningar 4.4 Insignalstabilitet . . 4.5 Stationära lösningar . 4.6 Resolventmatrisen . . 4. i Stationära lösningar och superposition . 4.8 Resonans . . . . . . . . . . . . . . 4.9 Instabilitet och linearitet . . . . . 4.10 Reella system av andra ordningen 4.11 Referenser . . . . . . . . . . . . .
61 61
5 Matrisfunktioner 5.1 Exponentialmatrisen . . . . . . . 5.2 Matrispotensserier . . . . . . . . 5.3 Lösning av tillståndsekvationer . 5.4 Matrispolynom och Cayley-Hamiltons sats 5.5 Bevis av stabilitetssatserna 5.6 Referenser . . . . .
85 85
102
6 Symmetriska matriser 6.1 Vektorgeometri . . . 6.2 Projektionsmatriser 6.3 Speglingsmatriser 6.4 Vridningsmatriser 6.5 Ortogonala matriser 6.6 Symmetriska matriser 6.7 Referenser . . .
105 105 107 108 109 110 112 118
7 Kvadratiska former 7.1 Matrisbeskrivning av kvadratiska funktioner 7.2 Taylorapproximation . . . . . . . . . 7.3 Diagonalisering av kvadratiska former 7.4 Klassifikation av kvadratiska former 7.5 Kvadratkomplettering . . . . . . . . . 7.6 Sylvesters tröghetslag och spektrumklyvning 7.7 Odämpade lineära svängningar 7. 7.1 Separation av variabler . . . . . . . . . 7. 7. 2 Diagonalisering . . . . . . . . . . . . . 7.8 Geometrisk tolkning av kvadratiska ekvationer
119 119 122 124
62 70
71 72 74
76 80 81 82 83
88 92 93
99
126 128 132 135 137 138 143
INNEHÅLL 7.9
Referenser
8 Numerisk lineär algebra 8.1 Preliminärt innehåll . 9
Insignal-utsignalmodeller 9.1 »Black box»-modeller av system . . . . . . . 9.2 Hur realistiska är insignal-utsignalmodeller? . Två typiska system . . . . . . . . 9.3 9.4 Linearitet och superposition .. . 9.5 Stegfunktioner och pulsfunktioner 9.6 Sampling . . . . . . . . . 9. 7 Kontinuerligt och diskret . . . . . 9.8 lmpulssvar . . . . . . . . . . . . . 9.9 Tidsinvarians och translationsinvarians 9.10 Stabilitet . . . . . . . . . . . 9.11 Kausalitet . . . . . . . . . . 9.12 Tillstånds- och JU-modeller . 9.13 Mer om systemmodeller . 9.14 Referenser
vii
147
149 149 151 151 155 156 158 159 162 163 164 169 173 177 178 182 182
10 Faltningar 10.1 Definition av faltning 10.2 Räkneregler för faltning . 10.3 Några standardfaltningar 10.4 Kausala faltningar . . . . 10.5 Faltningar och systemteori 10.6 Stegsvar . . . . . . . . . . 10. 7 Faltningar i sannolikhetsteorin 10.8 Referenser
185 185 188
11 Generaliserade funktioner 11.1 Varför räcker inte vanliga funktioner? 11.2 Deltafunktionen . . . . . . . 11.3 Allmänna distributioner . . . . . . . . 11.4 Derivation av distributioner . . . . . . 11.5 Likhet av distributioner. Integration .. 11.6 Konvergens. Fourierserier .. . 11.7 Faltning av distributioner .. . 11.8 Distributioner och systemteori 11.9 Passningsmetoden . . . . . . .
199 199 202 203 206 214 216 220 221 223
190
192 193 194 196 197
INNEHÅLL
viii
11.10 Referenser
226
12 Frekvensanalys 12 .1 Överföringsfunktionen . . . . . . . 12.2 Exempel på överföringsfunktioner 12.3 Frekvensfunktionen . . . . . . . . 12.4 Superposition och frekvenssyntes . 12.5 lmpulssvar och överföringsfunktion . 12.6 Frekvensanalys av ett värmeledningsproblem 12. 7 Olineära och tidsvariabla system . . . . . . 12.8 Poler, nollställen och överföringsfunktionen 12.9 Resonans . . . . . . .
229
13 Fouriertransformationen 13.1 Definition . . . . . . . 13.2 Grundläggande exempel . 13.3 Enkla räkneregler . 13.4 Derivationsreglerna . . . 13.5 Inversionsformeln . . . . 13.6 Bevis av Fouriers inversionsformel 13. 7 Faltningssatsen . . . . . . 13.8 Parsevals formel . . . . . . . . . . 13.9 Obestämdhetsrelationen . . . . . . 13.10 Fouriertransformation av deltafunktioner 13.11 Tvådimensionell Fouriertransformation . . 13.12 Hur Fouriertransformerar man i praktiken? 13.13 Stabila system i tids- och frekvensområdet 13.14 Referenser . . . . . . . . . . . . . . . . . .
253 253 254 256 257 260 265 269 271 273 276 282 283 285 287
14 Laplacetransformationen 14.1 Definition av Laplacetransformationen . . . . . . 14.2 Fourier- och Laplacetransformationerna . . . . . 14.3 Laplacetransformation och analytiska funktioner 14.4 Derivationsregler . . . . . . . . . . . . . . . . . . 14.5 Inversionsformeln . . . . . . . . . . . . . . . . . 14.6 Funktionsteoretisk härledning av inversionsformeln 14. 7 Laplacetransformation och distributioner . . . . . 14.8 Rationella transformer . . . . . . . . . . . . . . . . 14.8.1 lnverstransformation genom partialbråk.suppdelning 14.8.2 Inverstransformation med residykalkyl 14.9 Lösning av differentialekvationer . . . . . . . . . . . . . .
289 289 292 294 295 297 298 299 301 302 305 309
229 231 234
239 241 242 245 246 248
INNEHÅLL
14.10 14.11 14.12 14.13 14.14
ix
Lineära system och Laplacetransformationen Den ensidiga Laplacetransformationen . Begynnelsevärden och ensidiga problem Faltningsekvationer . . . . . . . . . Begynnelse- och slutvärdessatserna .
311
15 Stabilitetsteori 15.1 Stabilitetsbegreppet . . . . . . . . . . 15.2 Återkoppling . . . . . . . . . . . . . . 15.3 Metoder för stabilitetsundersökningar 15.4 Argumentvariationsprincipen . . . . . 15.5 Nyqvistdiagram . . . . . . . . . . . . 15.6 Rouches sats och algebrans fundamentalsats. 15.7 Poler på randen och i oändligheten . 15.8 Routh-Hurwitz kriterier 15. 9 Referenser . . . . . . . . .
325
A Skalära differentialekvationer
345
B Spektralsatsen B. l Bevis av spektralsatsen . . . . . . . . . . . . . B.2 Gram-Schmidts algoritm och QR-faktorisering
349 349 352
Referenser
355
Index
356
314 316
321 322
325 326 328 329 333
337 339 342
342
X
INNEHÅLL
1 Dynamiska system
1.1
Matematisk beskrivning av system
Ett system är en helhet, bestående av delar som ömsesidigt påverkar varandra. Denna kurs handlar om matematiska metoder för att beskriva och undersöka tekniska system. System 1 kan indelas i naturliga och av människan uppbyggda. Några exempel:
naturliga system
människoverk
• solen och de nio planeterna • population av rovdjur och population av bytesdjur
• radiomottagare • byggnad • Sveriges ekonomi
I stort sett kan man använda samma matematiska metoder för att studera båda dessa typer av system. Man kan aldrig fånga alla egenskaper hos verkligheten genom en matematisk (och inte heller genom någon annan sorts) beskrivning utan måste göra sig en förenklad bild, en modell, av de egenskaper som är relevanta. Vad som skall anses vara relevant bestäms av ändamålet med beskrivningen. Olika typer av prövning får sedan visa det praktiska värdet av denna modell. Vi skall här syssla med matematiska modeller av system. En sådan modell är uppbyggd av systemvariabler och av matematiska relationer mellan sådana variabler. Systemvariablerna, som vi kan beteckna med x 1 , x 2 , x 3 , ... mäter alla relevanta egenskaper hos systemet. I solsystemet är bland annat koor1
0m elektroteknikerna ger en sträng definition av krets,
1
så skall jag ge en av system.
2
KAPITEL 1. DYNAMISKA SYSTEM
dinater för himlakropparnas lägen och hastigheter vid en viss tidpunkt systemvariabler. i radiomottagaren har vi strömstyrkor och spänningar, temperaturer, lägen för de olika komponenterna. I regel befinner sig systemet i kontakt med en omvärld och man måste också ta hänsyn till dess växelverkan med denna. Denna växelverkan beskrivs med insignaler eller styrsignaler, som vi betecknar med w1 , w2 • w 3 ... samt med utsignaler eller mätsignaler, som betecknas med y1 • Y2, y3 •.... För radioapparaten kan t ex antennsignal och inställning av frekvens och volym betraktas som insignaler och ljudet från en högtalare som utsignal. En speciell typ av insignaler är mer eller mindre okontrollerbara störningar. Kopplingen mellan modell och verklighet är ett avgörande led men hör främst hemma i tillämpningsämnena. Man bör härvid tänka på att en framgångsrik modell i hög grad påverkar vår uppfattning om vad verkligheten är.
1.2
Statiska och dynamiska system
Ett viktigt led i systembeskrivningen är att bestämma hur systemvariablerna förändras med tiden. Denna förändring sker dels på grund av systemets inre egenskaper och dels på grund av påverkan från insignalerna. I detta avseende kan man skilja mellan två huvudtyper av system, nämligen statiska system och dynamiska system. I ett statiskt system beror systemvariablernas (eller åtminstone utsignalernas) värden vid en viss tidpunkt endast av insignalernas värden vid samma tidpunkt. Exempel på system som väl beskrivs av statiska modeller är rent resistiva elektriska nät och mekaniska system i jämvikt. Vanligare är att systemvariablerna vid en viss tidpunkt inte bara beror av insignalerna utan även av systemets tidigare historia. Ett sådant system kallas för ett dynamiskt system. En bil i rörelse måste beskrivas som ett dynamiskt system, ty dess hastighet och än mer dess läge beror (tyvärr) inte bara av styrsignalerna för ögonblicket utan också på utgångsläge och hastighet samt av styrningen i det förflutna. Elektriska nät som innehåller energilagrande komponenter, t ex spolar och kondensatorer, måste också beskrivas dynamiskt. Strängt taget är statiska modeller n~tan alltid approximationer av mer noggranna dynamiska, och används framför allt därför att de är mycket enklare att behandla matematiskt än de dynamiska. Skillnaden mellan dynamiska och statiska system återspeglas i de använda matematiska metoderna. För beskrivning av system, som definitionsmässigt består av flera komponenter, krävs alltid flera variabler. Härav följer att såväl flervariabelanalys som lineär algebra är nödvändiga hjälpmedel. För statiska system är detta också tillräckligt medan för dynamiska även tidsberoendet måste behärskas. Det vanligaste sättet att göra detta är med användning av ordinära differentialekvationer. För speciella system finns även andra metoder, tex integraltransformer och faltningar,
1.3. TILLSTÅNDSVARIABLER
3
som vi skall behandla i den senare delen av kursen. Utanför ramen av denna kurs faller flera viktiga typer av systemmodeller. Tidsdiskreta (samplade) system beskrivs med differensekvationer i stället för med differentialekvationer. Andra typer av system, från fältteori och kontinuerliga mediers mekanik, beskrivs med partiella differentialekvationer. I båda dessa fall är de metoder som behandlas i föreliggande kurs grundläggande, men måste kompletteras. En mera avvikande men i en hel del drag parallell teori förekommer i digitalteknik och datorteknik, där den använda matematiska teorin bygger på (ofta olineär) algebra. Se t ex Cedervall & Lindh ( 1995).
1.3
Tillståndsvariabler
Det finns i regel ett stort antal systemvariabler, men de är inte alla oberoende. Om man känner ett lämpligt antal av dem så kan värdena av de övriga beräknas. För att underlätta arbetet bör man välja ut en uppsättning systemvariabler som är så liten till antalet som möjligt men ändå ger all behövlig information. De som ingår i en sådan uppsättning kallas för tillståndsvariabler. Den grundläggande iden är att i tillståndsvariablemas värden vid en viss tidpunkt skall all den information om systemets tidigare historia vara sammanfattad, som behövs för att förutsäga systemets vidare utveckling, om insignalema är kända. Beskrivningen av system med hjälp av tillståndsvariabler har relativt nyligen vunnit en stor popularitet inom tekniska vetenskaper, bland annat därför att den är väl lämpad för behandling av system med hjälp av datorer. Den har dock en lång historia bakom sig och utgångspunkten finns i Lagranges och Hamiltons analytiska mekanik. En helt exakt definition blir med nödvändighet komplicerad, men följande beskrivning är tillräcklig för våra behov: Antag att en matematisk modell för ett system är given. En uppsättning (oberoende) systemvariabler kallas för tillstAndsvariabler för systemet om • värdena av samtliga systemvariabler vid en viss tidpunkt kan beräknas ur värdena av tillstdndsvariablerna och insignalerna vid samma tidpunkt. • systemets tillstånd vid varje tidpunkt t efter tiden t 0 är fullständigt bestämt av ( a) tillståndsvariablernas värden vid tiden t 0 och (b) insignalemas värden under tiden från t 0 till den aktuella tiden.
KAPITEL 1. DYNAMISKA SYSTEM
4
Som ett exempel på system kan vi ta en boll som kastas. Från eventuell rotatio~ rörelse (skruv) bortses i modellen. Styrsignalerna är de yttre krafterna, väsentligen tyngdkraft och luftmotstånd. Vilken information behövs för att vid ett visst ögonblick helt bestämma bollens framtida rörelse? Vi måste ha läget, beskrivet tex av lägeskoordinater x 1 , x 2 , x 3 , men detta är inte tillräckligt. Däremot är den framtida rörelsen entydigt fastlagd om vi förutom utgångsläget även känner utgångshastigheten. Som tillståndsvariabler kan vi alltså välja variablerna x 1 , x 2 , x 3 , x 4 , x 5 , x 6 , där x 1 , x 2 , x 3 är lägeskoordinater och x 4 = v1 , x 5 = v2 , x 6 = v3 är motsvarande hastighetskom ponenter. Det bör påpekas att valet av tillståndsvariabler alltid kan ske på många olika sätt. I exemplet med bollen kan man naturligtvis använda vilket koordinatsystem som helst i rummet för att beskriva läget, och även många andra modifikationer är möjliga. Vid konkreta systembeskrivningar är det sällan så enkelt att systemets utveckling direkt kan utläsas ur tillståndsvariablernas begynnelsevärden och insignalerna. I stället har man ett indirekt samband i form av ett antal ekvationer som knyter samman variablerna och ur vilka tillståndsvariablerna måste bestämmas. I föreliggande kurs rör det sig vanligen om ett antal differentialekvationer. I denna typ av systemmodeller förutsätts tillståndsvariablerna uppfylla ekvationer av formen dx1
dt
= /1(X1, X2, •. •, Xn, t)
dx2
= '2(x1, X2, ... , Xn, t)
dt
där / 1 , '2, ... , In är givna funktioner, bestämda av systemet och insignalerna. Ett system av differentialekvationer av en sådan form sägs vara i normalform eller tillståndsform. I teorin för ordinära differentialekvationer (som behandlas i kursen Olineära dynamiska system, set ex Spanne {1995b) och Andersson & Böiers {1992)) bevisas att under lämpliga villkor på funktionerna / 1 , '2, ... , /n har detta system alltid entydig lösning för givna begynnelsevärden. Härav följer att namnet tillståndsform är berättigat för sådana system av differentialekvationer.
1.4
Konstruktion av tillståndsmodeller
Konsten att ställa upp tillståndsmodeller för system, och matematiska modeller över huvud taget, är en del av tillämpningsämnena snarare än av matematiken.
TILLSTÅNDSMODELLER
5
Källorna till tillståndsekvationer är mångfaldiga. Man kan använda geometriska samband, allmänna fysikaliska principer som Newtons rörelseekvationer, Maxwells fältekvationer, balansekvationer som uttrycker massans och energins konstans men också speciella komponent- och materialekvationer som Ohms lag, Hookes lag etc. Vi skall här endast ge ett fåtal exempel för att konkretisera de ovan införda allmänna begreppen.
Exempel 1.1 (Resistiv spänningsdelare) I nätet i figuren väljer vi att betrakta spänning+ en w som insignal och spänningen y = v2 som utsignal. Betrakta systemvariablerna v 1 , v2 och w i. Dessa är inte oberoende utan Kirchhoffs spänningslag ger sambandet w = v1 + v2 och Ohms lag för resistorerna ger sambanden v 1 = R 1 i och v2 = R2 i. Härur fås direkt utsignalen som en funktion av insignalen, R2 y(t) = R R w(t). l
+
R
i
+
l
+
y
2
Utsignalens värde vid tiden t beror alltså endast av insignalens värde vid samma tidpunkt. Detsamma gäller för alla andra systemvariabler. och systemet är statiskt. För att beskriva det behövs överhuvudtaget inga tillståndsvariabler. En förutsättning för att modellen skall vara riktig är att spänningsdelaren (som i schemat) är obelastad, eller med andra ord att ingen ström går genom utgången. Denna förutsättning kom.mer vi att stillatigande göra i fortsättningen även i andra elektriska exempel. D
Vi skall nu se på vårt första dynamiska system.
Exempel 1.2 (Ett första ordningens lågpassfilter) Om resistorn R 2 byts ut mot en kondensator I med kapacitans C så ersätts Ohms lag med sam0 • + bandet Cd½ - . --1. dt w Vi skall se att i detta system kan vi välja v 2 som (ensam) tillståndsvariabel. Med x = v 2 får vi 0 W
eller
= Ri + V2 = RC: + X
R
+
V1
I C
_t+
T~
0
+ y
0
KAPITEL 1. DYNAMISKA SYSTEM
6
Sambandet mellan insignalen w och systemvariabeln x är nu en differentialekvation. Systemet är dynamiskt, eftersom dess historia inverkar på systemvariabler och utsignal. Fysikaliskt beror detta på att laddning kan lagras i kondensatorn. I detta enkla fall kan vi direkt skriva upp lösningen till differentialekvationen (se vidare appendix A). Den blir av formen
Om vi känner begynnelsevärdet x( t 0 ) och insignalen w( t) under tidsintervallet to ~ T ~ t så kan alltså värdet av x(t) bestämma.5. Detta bevisar att x är en tillståndsvaD riabel för systemet.
Exempel 1.3 (Ett hydrauliskt system) Systemet i denna figur påminner måhända om ett på LTH välbekant men förutsätts fungera. Kärlen I och Il är lika och har tvärsnittsarean A. Utifrån tillförs de vattenkvantiteterna w 1 (t) resp w 2 (t) per tidsenhet. Utströmningen ur kärlen förutsätts följa Torricellis lag; utströmningshastigheten är proportionell mot kvadratroten ur vattenhöjden i kärlen. Balansekvationerna volymändring ger sambanden
I
Xt
1~•V,•I-~----
inströmning
y
utströmning
~ = W1 ( t) - kv'Xi°
dV2 dt
och eftersom
= kv'Xi + W2(t) - ky'xi
Vi = x 1 A och V2 = x 2 A så får vi ett system av differentialekvationer,
TILLSTÅNDSMODELLER
7
Detta är ett system på tillståndsform där funktionerna f 1 och
I
h ges av
1
fi(x1, x2, t)
= A (-kv'XJ + w1(t))
h(x1, X2, t)
=:
(kv'XJ - k.jxi + w2(t))
och höjderna x 1 och x 2 kan användas som tillståndsvariabler. I detta exempel finns ingen enkel formel som direkt uttrycker tillståndsvariablema x 1 och x 2 med hjälp av begynnelsevärdena x 1 (t0 ) och x 2 (t 0 ) samt insignalerna w 1 (r) och w 2 (r), såsom fallet var i föregående exempel. Det visar sig att den enda möjligheten är en numerisk lösning av differentialekvationerna. Man skulle också kunna tänka sig att bygga en fysikalisk modell i mindre skala, men det är då osäkert i vad mån denna modell återspeglar det verkliga systemets egenskaper. (En sådan modell finns verkligen av det ovan åsyftade verkliga systemet. men den lär fungera □ i motsats till det stora systemet). I exemplen 1.1, 1.2 resp 1.3 ovan har vi använt 0, 1 resp 2 tillståndsvariabler. Antalet tillståndsvariabler är en viktig systemparameter, som inom vida gränser är oberoende av matematiska omformuleringar av systemmodellen. Detta antal brukar kallas för systemets ordning. Med ordningen av en systemmodell menas antalet tillståndsvariabler i modellen. Svårighetsgraden vid systemanalysen ökar nästan alltid snabbt med ordningen. I våra två första exempel lyckades vi finna ett explicit uttryck för utsignalen som funktion av insignal och begynnelsevärden. Detta berodde framför allt på att systemen var av första ordningen men också på att de var lineära. I systemet i exempel 1.3 finns ingen analytisk lösning, dvs ingen allmän formel för lösningen, utan numeriska metoder måste tillgripas. I praktiken uppkommer ofta system av mycket hög ordning. Ett elektriskt nät som innehåller p kondensatorer och q induktorer brukar vara av ordning p + q, och lämpliga tillståndsvariabler är ofta kondensatorspänningarna och induktorströmmarna. Ett mekaniskt system som innehåller N fria partiklar är av ordning 6N. Som tillståndsvariabler kan väljas 3N lägeskoordinater, 3 för varje partikel, och motsvarande 3N hastighetskomponenter. Solsystemet är alltså (om man försummar månar, små.planeter och kometer) av ordning 60 (solen + 9 planeter). En mol av en (enatomig) gas utgör i denna mening ett system av ordning ~ 6 • 6 • 1023 (6 • Avogadros tal). En direkt fullständig analys av ett system av så hög ordning är naturligtvis omöjlig, vare sig man räknar analytiskt eller använder dator.
8
KAPITEL 1. DYNAMISKA SYSTEM
Inom fysiken är det också vanligt med tillståndsmodeller som inte är av ändlig ordning. dvs som inte kan beskrivas med ändligt många tillståndsvariabler. Sådana modeller kommer från olika fysikaliska fältteorier och behandlas närmare i kursen Kontinuerliga system. Vid praktiskt arbete måste man alltid göra en avvägning mellan precisionen och enkelheten i de modeller man använder. En approximativ modell, kanske av lägre ordning. som är möjlig att analysera är i regel bättre än en mer exakt modell, där man inte kan lösa de ingående ekvationerna. Konsten att göra »ingenjörsmässiga approximationer» är mycket viktig. även om i dag ganska komplicerade modeller kan hanteras med numeriska metoder.
Inläsningsuppgifter: Kapitel 1 b) tillståndsvariabler 1 Redogör för lärobokens beskrivning av c) ordning av system. följande begrepp: a) matematisk modell
2 Lineära tillståndsmodeller
Detta kapitel innehåller en inledning till lineär modellering och lineära tillståndsmodeller. Detta ingår visserligen inte i kursen Lineära system, men det känns säkert mer meningfullt att lära in de matematiska metoderna om man sett några konkreta fall där de kan användas.
2.1
Lineära och olineära funktioner
Vid matematisk modellering av system uppkommer funktionssamband mellan flera variabler av typen Yl= /1(X1,X2, ... ,Xn) Y2 = h(x1, X2, ... , Xn)
eller i vektorform y
= f(x)
Allmänna sådana samband kan vara mycket komplicerade. I praktiken finns det olika sätt att beskriva dem, beroende på om de tagits fram experimentellt, genom numeriska beräkningar eller genom helt teoretiska överläggningar: • genom analytiska formler • grafiskt, med kurvor • med tabeller 9
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
10
Exempel 2.1 Sambandet
är av den form som brukar förekomma i övningsexempel i flervariabelanalysen. Det är klart olineärt. D Ett exempel av mer direkt praktisk betydelse är följande.
---------
Exempel 2.2 (Transistor) En transistor (it ex GE-koppling) kan vid måttligt höga frekvenser representeras som en resistiv fyrpol, där endast två av de fyra variablerna (i 1 , vi, i 2 , v2 ) är oberoende. Vanligt är att välja i 1 och v2 som oberoende variabler. vilket ger funktionssamband av formen
1
i2
+ V1 I •
- - - - - - - - ... Z2
(hybridmodell av transistor) där h 1 och h2 är olineära funktioner. Det finns inga enkla analytiska uttryck för dessa funktioner.
-8 -100 A i1
i1
=Is= konstant
=Is= OµA
0
-10
(V)
I datablad för transistorer representeras de oftast genom kurvor som t ex ger i 2 som funktion av v2 vid olika konstanta värden på i 1 , alltså en grafisk representation av h2, och motsvarande kurvor för h 1 . D
+
2.1. LINEÅRA OCH OLINEÅRA FUNKTIONER
11
Grafisk representation är mycket använd då man har sammanlagt högst tre eller fyra variabler men svårigheterna med en sådan blir stora så snart man har med mer än två oberoende variabler att göra. Instängda som vi är i ett tredimensionellt rum har vi ingen större förmåga att grafiskt föreställa oss samband mellan mer än tre variabler. Också datorer har än så länge svårigheter att använda data givna i form av diagram. Det kan t ex nämnas i detta sammanhang att mycket arbete läggs för närvarande ner (även vid LTH) på att finna enkla analytiska olineära transistormodeller att använda vid datorsimulering av integrerade kretsar. En tabellrepresentation blir likaledes ogörlig för funktioner av många variabler. Antag att vi vill beskriva en funktion y = f(x 1 , x 2 , x 3 ) av tre variabler i intervallet -10 ~ xi ~ 10. i = 1, 2, 3 med varje Xi givet med 2 decimaler. Då krävs ungefär 20/0.01 = 2 -1Q3 värden av varje Xi, dvs totalt (2 • 1Q3) 3 ~ 1010 funktionsvärden. På en tabellsida får det plats högst ca 1000 värden, och tabellen skulle uppta ungefär 10000 volymer med vardera 1000 sidor. Detta verk skulle ta upp större plats än matematiska institutionens hela bokbestånd. Situationen blir en helt annan om man inskränker sig till lineära samband, alltså samband av formen:
=
Yl Y2
=
a11X1 a21X1
+ +
a12X2 a22X2
+ ··· + + ··· +
a1nXn a2nXn
eller i matrisform y=Ax där A är matrisen
A=
au
a12
a1n
a21
a22
a2n
llm1
llm2
llmn
(Se Sparr (1995).) Lineära funktioner är fullständigt bestämda av sina matriser A. Ett lineärt samband mellan n oberoende variabler x; och m beroende variabler Yi kräver alltså lagring av endast mn tal ai;, vilket ju är ett måttligt antal även för ganska stora värden på m och n. Det är inte bara utrymmesproblem som gör att man gärna arbetar med lineära modeller. För lineära problem kan man använda kraftfulla teoretiska resultat ur den lineära algebran. Det finns numera också mycket effektiva metoder att till måttlig kostnad lösa lineära problem på dator. Även små olineära problem leder lätt in i ett matematiskt okänt (men mycket intressant) område.
12
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
2.2
Källor till lineära modeller
I denna kurs skall vi nästan uteslutande behandla matematiska metoder för lineära system. Källorna till lineära systemmodeller är väsentligen två: • lineära grundekvationer • linearisering av olineära samband Många fysikaliska lagar är från början lineära. Ett exempel är Kirchhoffs lagar i strömkretsteorin, som vi behandlar i avsnitt 5. Strömkretsteorin ger upphov till helt lineära ekvationer om ingående komponenter tillräckligt väl kan beskrivas med lineära ström-spänningssamband. Detta är ofta fallet med resistorer, kondensatorer och induktorer. För dessa brukar man använda följande standardmodeller:
R
C
L
-'VV'v-
~~
~
V= Ri
Cdv dt
=i
di v=Ldt
En diod är däremot en typiskt olineär komponent. Många dioder kan modelleras med ström-spänningssambandet z
i= io (e""i.UT - 1) V
där i 0 , Å, q och k är fysikaliska konstanter. En transistor kan uppfattas som sammansatt av dioder, och den ovan antydda transistormodellen är givetvis olineär.
2.3
Linearisering
Om ett samband från början är olineärt, så kan det under vissa förutsättningar approximeras med ett lineärt. Hela differentialkalkylen bygger på denna ide, nämligen att en (deriverbar) kurva nära en punkt kan ersättas med sin lineära tangent. Denna procedur kallas för linearisering.
2.3. LINEARlSERING
13
Antag att vi har ett samband
Yl= Y2
=
/1(X1,X2, ... ,Xn)
h(x1, X2,
... , Xn)
eller i vektorform y
= f(x)
som är differentierbart i punkten x" = (xt, x~, ... , x:). Inför avvikelserna 6x och 6y från vilopunkten eller arbetspunkten (x", y 0 ) = (x", f (x")) genom {
6x 6y
= x- x° = y - y = f(x) 0
f(x°)
Definitionen av differentierbarhet säger då att sambandet mellan 6x och 6y är lineärt på ett litet fel när: 6y ~ A6x om 6x är tillräckligt litet. Matrisen A kallas för funktionalmat.risen för f och erhålles som 8/i A = f'(x
0 )
8x1
8/i 8Xn
8/m 8x1
8/m 8xn
=
där derivatorna beräknas i punkten x = x". Mer precist är feltermen 6y - A6x liten i jämförelse med 16xl då 6x är litet, det vill säga l6y l~x16xl
➔0
då 6x
➔0
Detta brukar uttryckas så att 6y
= A6x + o(l6xl) då 6x ➔ 0.
För denna beteckning (litet ordo) hänvisas till Persson & Böiers (1988), sid 374-377. Om avvikelsen 6x från arbetspunkten är tillräckligt liten så är det rimligt att strunta i resttermen. Vi får då ett lineärt samband mellan de nya variablerna 6x och 6y: 6y = A6x
KAPITEL 2. LINEÄRA TILLSTÅNDSMODELLER
14
Denna relation kallas för lineariseringen av det olineära sambandet y = f (x) kring arbetspunkten x = x 0 • Proceduren beskrivs i tillämpningslitteratur ofta så att man »serieutvecklar / och stryker högre ordningens termer». Detta är matematiskt sjuttonhundratalsspråk, och linearisering har mycket lite att göra med serieutveckling. För att få fram koefficienterna i funktionalmatrisen kan man i konkreta fall förfara på en rad olika sätt: • analytisk beräkning av derivator • numerisk beräkning av derivator • mätningar i kurvor (derivata
= lutning)
• mätningar direkt i det fysikaliska systemet Vi ger först ett exempel på analytisk linearisering.
Exempel 2.3 Linearisera funktionen {
kring arbetspunkten x 0 Lösning: Eftersom
= ez1z2 2 2 Y2 = X1 + X2 Yl
= (1, 2).
så blir funktionalmatrisen
och A = f'(xo)
= [2e 2 2
Alltså är
6y
= [ 1~:~
e2 ] 4
= [ 14.78
!:: ]
7.35 2.00 4.00
/;x
l
+ o(l6xl)
Vi jämför det exakta och det lineariserade värdet på 6y för några värden på 6x:
15
2.3. LINEARlSERING
6x1 0.01 -0.08 -0.10 -0.10 0.20 2.00 3.00
6x2 0.00 0.04 0.04 0.10 0.10 1.00 5.00
ÖYtez
0.15 -0.86 -1.86 -0.77 5.04 8096 1.4. 1012
6yuin
Ö'Y'iez
ÖY2lin
0.15 -0.898 -2.22 -0.74 3.69 37 81
0.02 0.01 -0.58 0.22 0.85 13 60
0.02 0.00 -0.60 0.20 0.80 8 26 D
Allmänt kan lineariseringsfelet o(l6xl) enligt Taylors formel uppskattas med Cl6xj 2 om f är två gånger kontinuerligt deriverbar. Konstanten C beror av andraderivatornas storlek för betraktade 6x. I exemplet ovan är konstanten C av storleksordningen 100 för små x. Ett typiskt exempel på linearisering hämtar vi ur ett kompendium i Tillämpad elektronik.
Exempel 2.4 {Småsignalanalys av en transistor) Då man använder en transistor för att förstärka analoga signaler vill man i regel att utsignalen skall ha samma form som insignalen, bara multiplicerad med en (stor) konstant. Detta åstadkoms genom att man anlägger förspänningar så att transistorn arbetar i en viss vilopunkt och sedan låter de önskade signalerna representeras av små avvikelser från denna vilopunkt. Den olineära transistormodellen från exempel 2.2 har formen
För små. avvikelser från vilopunkten
(= (isQ,VBEQ,icQ,VcEQ) i elektroniken, där Q står för quiescent de lineariserade sambanden
= vilande) får vi
KAPITEL 2. LINEÄRA TILLSTÅNDSMODELLER
16
Elementen i funktionalmatrisen kallas för h-parametrar och de vanliga beteckningarna är
[~::
~~] = [~;
~:]
De kan erhållas ur datablad. grafiskt eller i tabellform. Exempelvis är ho
= h 22 = 8h2 = ( 8i2 ) 8v2
8v2
11
lika med lutningen av v2 - i 2-kurvan (i 1 = konstant) i arbetspunkten. De komplicerade olineära sambanden ersättas alltså med de enkla lineära
{
6~1
= hj6i~ + hr6V2
6z2
= h16z1 +
h0 6v2
som kan användas då transistorn befinner sig tillräckligt nära arbetspunkten.
D
Ett annat viktigt exempel på linearisering är analysen av ett mekaniskt system som utför små svängningar kring ett arbetsläge. Vi kommer att behandla detta fall utförligare längre fram.
2.4
Linearisering av tillständsmodeller
Ett system av tillståndsekvationer dx dt
sägs ha en jämviktspunkt i x att den konstanta funktionen x
= f(x. t)
om f(i, t) = 0 för alla t. Detta innebär precis = i är en lösning till systemet.
=i
Exempel 2.5 Fontänen i exempel 1.3 beskrivs av funktionerna /1(X1,X2,t) { h(x1, x2, t)
! =! =
(-kv'xI +w1(t)) (kv1Xi - k./x2 + w2(t)).
Antag att insignalerna (tillrinning) w 1 och w 2 är konstanta i tiden. Ett jämviktstillstånd har vi precis då vätskenivåerna är konstanta, och alltså inflöde och utflöde precis balanserar varandra hela tiden. Vi skall bestämma möjliga jämviktspunkter (x 1, x2) = (i: 1, i: 2). De ges av jämviktsekvationerna
2.4. LINEARISERING
17
Alltså finns precis en jämviktspunkt, nämligen (i- 1 , i 2 )
= (( wk1 )2' (w1 +k W2) 2)
Om man startar med precis dessa vattennivåer, så kommer inrinning och utrinning att exakt balansera varandra. Om man startar med andra nivåer, så kommer vätskeytornas höjd att variera med tiden. D Man arbetar ofta nära ett jämviktsläge. Det är då naturligt att införa avvikelserna från detta som nya tillståndsvariabler och linearisera funktionerna f kring jämviktsläget. Med får vi f(k
+ 6x, t)
~ f(k, t) +f'(k, t)6x
....__,_, =O
och tillståndsekvationerna kan approximeras med lineära ekvationer för 6x,
d~;) =
f'(k, t)6x.
Exempel 2.6 För fontänen får vi genom derivation funktionalmatrisen
f' (x) = [ :::
::: ] =
8h 812
8x1
~ ~ 2A
[-
_1_
vXi
8x2
och insättning av värdena i jämviktspunkten ger
f'(k)
= !i:._ [- :11 2A W1
Detta ger det lineariserade systemet
I
d( 6x1) dt
= _ }!:__ _!__ 6x1
d( 6x2)
= }!:__ _!__6x1 -
dt
2Aw1
2A W1
vilket också kan skrivas med matriser som
}!:__ 2A W1
1
+ W2
6x2
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
18
Här beror alltså tidsderivatorna av tillståndsvariablerna lineärt av variablerna själva. I nästa kapitel skall vi se att varje system av denna typ har en relativt enkel D explicit lösning. i motsats till det ursprungliga olineära systemet. Vad är nu värdet med de lineariserade systemen? Det är ju inte självklart att ekvationer som liknar varandra har lösningar som liknar varandra, och det finns många praktiska exempel på motsatsen. För system av tillståndsekvationer finns det dock en hel del resultat som säger att lösningarna till det olineära systemet och lösningarna till det lineariserade systemet ligger nära varandra, om de inte avviker för mycket från jämviktspunkten. Detta är dock en fråga som vi inte kan avgöra här utan som hör hemma i kursen i olineära dynamiska system. Set ex Spanne (1995b) och för en mer avancerad framställning Guckenheimer & Holmes (1983).
2.5
Lineära tillståndsmodeller
Vi skall en längre tid framåt endast arbeta med systemmodeller av en speciell typ, nämligen lineära, tidsinvarianta tillståndsmodeller av ändlig ordning. Dessa modeller har definitionsmässigt följande matematiska struktur: Systemet har • n tillståndsvariabler x 1 , x 2 , · • p insignaler w1 , w2 , · ·
· ,
· · • Xn
w,,
• q utsignaler Y1, 112, · · · , Y9
Sambandet mellan tillstånd och insignaler ges av ett system av differentialekvationer dx1
dt
= a11X1 + a12X2 + · · · + a1nXn +
dx2
= a21X1 + a22X2 + · · · + a2nXn + ~1 W21 + · · · + ~,,w,,
dt
dxn dt
= On1X1 + On2X2 + · · · + OnnXn +
bu W1
bnl W1
+ · · · + b1,,w,,
+ · · · + bnpWp
och utsignalerna bestäms av tillståndsvariabler och insignaler genom ett statiskt lineärt system Y1 Y2
= C11X1 + C12X2 + · · · + C1nXn + d11w1 + · · · + d1pWp = C21X1 + C22X2 + · · · + C2nXn + d21 W1 + · · · + d2pwp
2.5. LINEÅRA TILLSTÅNDSMODELLER
19
Parametrarna Cli;, bi;, Ci; och ~; antas vara oberoende av tiden t (tidsinvarians). I matrisform får systemekvationerna det betydligt enklare utseendet dx
{
d = : =
Ax+Bw Cx+Dw
där A, B, C och D är konstanta matriser. Den viktigaste av dessa är A. som ofta kallas för systemmatrisen. Den beskriver systemets inre dynamik, hur de olika tillståndsvariablerna är dynamiskt kopplade till varandra. Lineära systemekvationer av typen ovan används för att beskriva ett mycket stort antal tekniska system. Anledningen till detta är framför allt att det finns effektiva metoder att lösa sådana ekvationer, alltså att bestämma hur tillståndet x varierar med tiden, då insignalen w är given. Vi skall ägna de närmaste kapitlen åt att ta fram sådana metoder. En annan anledning till att modellerna ovan är så användbara är att man genom att välja en lämplig ordning n och sedan anpassa A, B, C och D till det aktuella systemet kan någorlunda tillfredställande simulera även system som egentligen har en helt annan struktur. Den verkliga världen är inte lineär och tidsinvariant men det är ofta möjligt att styra den om man låtsas det. Vi övergår nu till två relativt enkla men typiska exempel, ett för fysiker och ett för elektriker.
Exempel 2.7 (Endimensionell diffusion) Neutronerna i ett stycke uran rör sig slumpmässigt i alla riktningar, vilket medför en nettotransport från områden med högre koncentration till områden med lägre koncentration. Vi skall ge en enkel modell för hur denna transport genom diffusion på.verkar koncentrationens tidsberoende. Betrakta en stav av längd l. Antag att neutronerna bara rör sig i stavens längsriktning. För att beskriva koncentrationens rumsberoende delar vi först staven i två delar, av längd h = l/2 och kallar genomsnittskoncentrationen i dessa för u 1 respektive u2 . Den mängd neutroner Am som från en del passerar en gränsyta under en viss tid är approximativt proportionell mot koncentrationen u i denna del och mot tidsintervallets längd samt omvänt proportionell mot delens längd h ( eller mot den sträcka som neutronerna i genomsnitt har att genomlöpa innan de åker ut), alltså:
f)
A
Am~DhuAt
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
20
där Där en proportionalitetskonstant. diffusionskonstanten, och A betecknar tvärsnittsarean. Denna regel är en variant av Ficks lag och kan härledas med statistiska resonemang. Eftersom u = m/Ah (total mängd 1 volym) så blir motsvarande koncentrationsändring D au= --uat h2 Utanför staven antas koncentrationen vara u = 0. !\ettoändringen av u 1 sammansätts då av tre bidrag
-0;] vilket ger differenskvoten
au 1
at : : :
D h2 (-2u1
+ U2)
Samma resonemang leder till en ekvation för u 2 , och en tänkt gränsövergång at --+ 0 leder till systemet
eller i matrisform
du dt
= Au,
A
i i ]=
= [- 2 D
_2D h2
h2
D [ -2 1 h2 1 -2
l
Detta är ett lineärt system av andra ordningen på tillståndsform, med u 1 och u 2 som tillståndsvariabler. Hittills har vi förutsatt att inga neutroner tillförs utifrån, nybildas eller sönderfaller. En nettonybildning w per tids- och volymsenhet ger upphov till två nya termer i systemet. Detta övergår i
d:ii { du
2
dt
=
-2i
=
D h2
U1
+
i
U2
+ W1(t)
U2
+ W2(t)
D
U1 -
2 h2
Detta är ett system på tillståndsform med ,.insignal,. w. Vi kan skriva det i matrisform
du dt
= Au+ Bw
2.5. LINEÅRA TILLSTÅNDSMODELLER
21
om vi tar 8 som enhetsmatrisen /. Som observerbar utsignal kan vi till exempel se på det totala utflödet per tidsenhet ur staven, vilket ges av
eller i matrisfonn y= Cu+Dw
där
o = [o o]. Modellen där vi delat in staven i endast två delar är naturligtvis mycket grov, men vi skall så småningom se att den ändå återspeglar viktiga egenskaper hos det fysikaliska systemet. En bättre modell kan man få genom att dela in i fler, n stycken lika delar, som då vardera har längden h = l / n. Exakt samma resonemang som ovan ger systemet
dun dt där d
=
i,
eller
du dt
= Au+ w,
-2 1 0 1 -2 1
0
0 0
0 0
0 0
0 1 -2
Detta är ett system av ordning n med tillståndsvariablema u 1 , u 2 , ... , Un. Av detta exempel framgår att samma fysikaliska system kan modelleras med en mängd av olika tillståndsmodeller, här av godtyckligt hög ordning. I förbigående kan vi anknyta till den kommande kursen Kontinuerliga system. Det är naturligt att förfina indelningen och låta n ➔ +oc. Då blir u en funktion av
22
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
en rumsvariabel x längs staven. u = u(x, t). Differenskvoterna övergår i rumsderivator: 1 ~u h 2 (u1r-i - 2u1r + u1r~d ➔ 2
ox
Systemet beskrivs nu av en partiell differentialekvation
au
~u -=D-+u·
(Soc:)
lJt
{)x2
med randvillkor u(O. t) = u(l, t) = 0. Den matematiska behandlingen av den partiella differentialekvationen uppvisar många analogier med den vi kommer att ge för system av ordinära differentialekvationer, men den kommer att behandlas först i Kontinuerliga system. Det bör dock nämnas att vid numerisk lösning av (S00 ) så går man motsatt väg mot den vi gått och approximerar derivatorna med differenskvoter. det vill säga ersätter (Soc:) med ( Sn) för ett tillräckligt stort n. Detta sker med finita differensmetoder och finita elementmetoder. D Vi skall nu försöka tillfredställa de läsare. som är mer intresserade av rent elektriska problem.
Exempel 2.8 (Wienfilter) I praktiska oscillatorkopplingar ingår ofta nätet i figuren, Wienfiltret. Vi väljer att använda kondensatorspänningarna v1 och v2 som tillståndsvariabler. Ur Kirchhoffs lagar erhålls sambanden +
+
w
+ Enligt definitionen av komponenterna gäller grenekvationerna
Ur dessa samband eliminerar vi de övriga systemvariablerna i, vRi, i R 2 och ic2 • Detta ger efter en stunds räkning {6 ekvationer - 4 eliminerade variabler= 2 ekvationer):
-dv 1 = - - -1V 1 dt
R1C1
dv,z dt
= - R1C2 Vi -
1
1
1
--V2
R1C1 (
1
R1C2
+ --W R1C1 1
+ R2C2
)
v2 +
1
R1C2 W
2.5. LINEÅRA TILLSTÅNDSMODELLER
23
vilket är ett system på tillståndsfonn med två tillståndsvariabler v1 och v2 och en insignal w. De i systemet ingående matriserna är
A= Antag nu att man av någon anledning är speciellt intresserad av en viss systemvariabel, till exempel spänningen över motståndet R 1 . Vi kan då betrakta denna som utsignal y och får sambandet
Även detta är lineärt med matriserna C=
[-1 -1],
= [1].
D
Totalt har vi alltså fått ett system av typen
{
!;=Av+ Bw y
= Cv+ Dw
med de ovan angivna konstanta matriserna A, B, C och D.
D
Lägg märke till att systemen i de båda föregående exemplen har exakt samma matematiska struktur, men att den fysikaliska tolkningen av ingående variabler och koefficienter är helt olika. I de närmast följande kapitlen skall vi ta fram ett antal metoder att lösa system med en sådan matematisk struktur. Lyckligtvis har sedan samma ekvationer alltid samma lösningar, oberoende av vilken den konkreta tolkningen är. Då vi tar fram den matematiska tolkningen är det i många fall onödigt att arbeta med den fullständiga tillståndsformen ovan. Till en början behöver man inte ange vilken utsignalen y är. Den fås ju omedelbart om insignalen och tillståndet är kända. Sedan spelar den speciella formen av koppling mellan insignal och tillstånd, vilken ges av matrisen 8, ingen större roll vid lösningen. Vi kan därför formulera det allmänna matematiska problemet på följande sätt: Vi vill lösa systemet
I~
=Ax+ f(t)
I
av differentialekvationer, där f(t) är en given vektorfunktion och A en given matris. Vi kommer i fortsättningen att kalla även ett sådant system för ett lineärt system av ( tidsinwrianta) tillstAndsekvationer.
24
2.6
KAPITEL 2. LINEÅRA TILLSTÅNDSMODELLER
Ekvationer med derivator av högre ordning
Vid modellering av tekniska system uppkommer ibland tidsderivator av högre ordning. Newtons rörelselag förknippar t ex kraft och acceleration, och innehåller alltså en tidsderivata av andra ordningen. Sådana system kan med en enkel omskrivning överföras på ekvationer av tillståndsform. Denna omskrivning används nästan alltid då man vill lösa ekvationerna numeriskt, och likaså då man vill visa existens- och entydighetssatser för deras lösningar. Vi beskriver metoden i ett exempel.
Exempel 2.9 Vi skall överföra differentialekvationen d3y dt3
.) = det(>./ -
8) = det(>./ -
s- 1AS)
s- 1AS. Då är
=
= det(>.S- 1/S - s- 1AS) = det(S- 1(>./ - A)S) = = det 5- 1 det(>./ - A) det 5 = det(>./ - A) = PA(>.)
ty detS- 1 detS = det(S- 15) =det/= 1.
D
Det karakteristiska polynomet är alltså bestämt av systemet självt och beror inte av den speciella uppsättning tillståndsvariabler som vi valt att använda i vår modell. Speciellt är egenvärdena entydigt bestämda av systemet. Detta är väl förenligt med den tolkning av egenvärden som svängningsfrekvenser och dämpningar, som vi gjort ovan. Det karakteristiska polynomet kan faktoriseras som
PA().) = (). - >.1)(>. - >.2) · · · · · (>. - >.n) = = )." - (>.1 + >.2 + · · · + >.n)>.n-l + · · · + {-l)">.1>.2 · · · · · ).n· Koefficienterna i PA kan direkt uttryckas i egenvärdena med hjälp av de kända sambanden mellan rötter och koefficienter i en polynomekvation. Två av dessa k. - au PA(>.) = det(>./ - A) =
-a21
-a12 >. - a22
-0n1
-an2
>. - Onn
= >." - (au+ a22 + · · · + Onn)>."- 1 + · · · + {-1)" det A. Den ena av dessa koefficienter är lika med determinanten och den andra är lika med summan av diagonalelementen. Denna senare har ett eget namn, spåret av A, och betecknas ofta med tr A (av engelska trnce = spår): DEFINITION
3.5 Med spåret tr A av en matris menas summan av diagonalelemen-
ten: tr A = au
+ a22 + · · · + llnn•
Vi ser att två av koefficienterna i det karakteristiska polynomet direkt kan uttryckas med hjälp av spår och determinant:
PA(>.) = >." - tr A>."- 1 + · · · + (-1)" det A. Speciellt följer att likformiga matriser har samma spår och samma determinant. Jämförelse mellan två olika uttryck för PA ger de viktiga sambanden
3.10. FOURIERS METOD
57
Sats 3.13 Om matrisen A har egenvärdena
Å1 , ... , Ån
så är
n
tr A =
Å1
+ Å2 + · · · + Ån
=
L
Åj
j=l
n
detA
= Å1Å2 ·····Ån= Il Å;. j=l
Dessa formler är av stor teoretisk betydelse och kan också med fördel användas vid kontroll efter beräkning av egenvärdena till en matris.
Exempel 3.11 Bestäm summan av egenvärdena till matrisen 13 97 85 ] A = [ 6 -2 21 44
Lösning:
Å1
+ Å2 + Å3 =
trA
=
13 + (-2)
6
1r
+6 =
17
D
Två matriser som är likformiga har samma karakteristiska polynom, men omvänd-
ningen behöver inte vara sann. Exempelvis har
båda det karakteristiska polynomet PA(Å) = Ps(Å) = Å2 , men 8 är inte diagonaliserbar Uämför exempel 3.8) och alltså ej likformig med diagonalmatrisen A. Icke diagonaliserbara matriser utgör dock det enda undantagsfallet. Om A och 8 är diagonaliserbara och PA = p8 så har de samma egenvärden och är alltså båda likformiga med samma diagonalmatris, därmed också med varandra.
3.10
Fouriers metod
För att anknyta framåt till kursen Kontinuerliga system skall vi här ge en alternativ besk.rivning av diagonaliseringsmetoden i avsnitt 3.6. Vi vill som vanligt lösa ett begynnelsevärdesproblem
dx
dt =Ax+ f(t),
x(to)
= Xo-
Om vi känner en bas Si, s-i, ... ,Sn som består av egenvektorer till A, så kan vi förfara på följande sätt.
58 KAPITEL 3. DIFFERENTIALEKVATIONER OCH SPEKTRALTEORI
Ctveckla den sökta lösningen x(t) och de givna data Xo och f(t) efter basen St, Si . ... ,Sn. dvs ansätt n
= I: x,.(t)s,.
x(t)
k=l
och bestäm Xok och
J,.
så att
n
f..,. är egenvektorer)
k=l
Insättning i differentialekvationen och begynnelsevillkoret ger n d' n
L ;,. s,. = L(>..,.x,. + j,.)s,. k=l
k=l
n
n
k=l
k=l
L x,.(to)s,. = L
Xo1.S1.-
Eftersom vektorerna s1 , Si, ... , 5n utgör en bas och alltså är lineärt oberoende, så måste koefficienterna på bägge sidor vara lika. Detta innebär att
dx,. dt
.
= >..,.x,. + /1.,
i1.(to)
= iot-
Detta är exakt samma ekvationer som vi tidigare härlett med matrisräk.ning och lösningen blir densamma. I teorin för kontinuerliga system måste man ofta arbeta med system som inte är av ändlig ordning. I stället för systemmatrisen A uppträder vanligtvis en partiell differentialoperator. Om man kan finna ett fullständigt system av egenvektorer (i differentialoperatorfallet kallade egenfunktioner) så går metoden ovan igenom på formellt exakt samma sätt. Den enda skillnaden är att den övre summationsgränsen blir oändlig, och den ändliga summa vi har i föreliggande kurs ersätts med en oändlig serie. Härvid uppkommer naturligtvis problem med konvergens, men dessa lämnar vi med glädje till senare kurser.
3.11. REFERENSER
3.11
59
Referenser
1. Åström (1976): kursbok i reglerteknik.
2. Spanne (1995a): lärobok i fortsättningskursen i matristeori. In~ehålbler bl a en beskrivning av vad man kan göra då matrisen inte är diagonahser ar. 3. Strang (1988): Lärobok i lineär algebra på ungefär_ samma ni~ som denna kurs, skriven av en av mästarna i tillämpad matematik och numenska metoder och dessa områdens pedagogik. 4. Strang (1986): En mycket personlig introduktion till tillämpad matematik (mest lineär).
Inläsningsuppgifter: Kapitel 3
1 Redogör för definitionen av a) egenvärde och egenvektor b) karakteristiskt polynom c) spektrum d) likformighet e) diagonalisering f) enkelt och multipelt egenvärde g) spår av matriser. 2 Vilka typer av matriser kan ha egenvärden och egenvektorer? 3 Vilket samband finns mellan egenvärdena och det karakteristiska polynomet för en matris? Härled detta samband. 4 Ange hur man kan diagonalisera en matris med hjälp av dess egenvektorer. 5 Ge ett exempel på en matris som inte är diagonaliserbar. 8 Finns det en icke diagonaliserbar matris
med enbart enkla egenvärden. 7 Finns det en diagonaliserbar matris med multipla egenvärden? 8 Ange och härled sambanden mellan spår, determinant och egenvärden för en matris. 9 Redogör för superpositionsprincipen (J) för lösningar till homogena lineära ekvati t 0 , beror både på Xo och på f(t), t > t 0 . Instabilitet innebär att måttliga förändringar av Xo eller av f leder till stora ändringar i x. Vi skall först studera fallet då insignalen är identiskt noll, f(t) = 0, t > t 0 • Systemet rör sig då fritt ( autonomt = självstyrande), och man brukar tala om fria svängningar (även då rörelsen egentligen inte är en svängningsrörelse). Vi gör följande definition:
4.2. FRIA SVÅNGNINGAR
63
DEFINITION 4.1 Systemet
dx dt
(H)
= Ax
'
A konstant
sägs vara
• stabilt om varje lösning till (H) har gränsvärdet O då t
➔
+oo.
• neutralt stabilt om alla lösningar är begränsade för stora t, men det finns lösningar som inte har gränsvärdet O.
• instabilt om det finns lösningar som är obegränsade för stora t. Matrisen A sägs vara stabil etc om systemet (H) har motsvarande egenskap.
Den ovan definierade typen av stabilitet skulle man kunna kalla för begynnelsevärdesstabilitet, ty den innebär att om man rubbar systemet från jämviktsläget x = 0, alltså ger det ett begynnelsevärde Xo i= 0, så återgår det fritt mot jämviktsläget. Hur kan man avgöra om ett system av typen (H) är stabilt? Vi måste se på den allmänna lösningen och startar med fallet då A är diagonaliserbar. Den allmänna lösningen är då som bekant av formen
där Å1, Å2, ... , Ån är egenvärdena till systemmatrisen och s1, s-i, ... , 5n är en motsvarande bas av egenvektorer. Eftersom le'tl = e 0
Bevis: Om ett egenvärde Å1c har positiv realdel. så är den för stora t obegränsade funktionen x(t) = eÅ"'s1c en lösning, och systemet är alltså instabilt. Om egenvärdet Å1c har realdelen 0, så är termen x(t) = eÅ"'s1c begränsad men har ej gränsvärdet 0. Om slutligen Re Å1c < 0, så är har denna term gränsvärdet O. Detta täcker nu alla fall i satsen. □ Stabiliteten av ett system beror alltså i det diagonaliserbara fallet enbart av var systemmatrisens egenvärden ligger i det komplexa (frekvens-)planet. Om alla egenvärden ligger strikt i vänster halvplan, så är systemet stabilt, annars inte. Typiska egenvärdesdiagram i de olika fallen är w
w
s
X
X X
X
s X
X
O"
O"
X X
w
s
O"
X X
X
X
stabilt: neutralt stabilt: instabilt: alla egenvärden i vänster alla egenvärden i Re s :$ något egenvärde i höger halvplan Re s < 0. 0, några på imaginära halvplan Res > 0. axeln. Vi ger nu exempel på system av olika typer. Först kommer ett stabilt system med reella egenvärden.
4.2. FRIA SVÄNGNINGAR
65
Exempel 4.1 Om w
A=
[ -21
1 -2
l
-3
-1
så är Å1 = -3, Å2 = -1 (se exempel 3.1) och alltså u(A) = -1. Den allmänna lösningen till (H) är (t) =c1e -t +c2 e -3t { x1 -t -3t X2 (t) = C1e - C2e och den har gränsvärdet O då t --+ +oo. Systemet är stabilt. Alla lösningar är lineärkombinationer av reella exponentialfunktioner, och inga egentliga svängningar kan förekomma. □ Nästa exempel är ett stabilt system med komplexa egenvärden.
Exempel 4.2 Om w
-2+ i X
A-[ 0 1] -5 -4
O"
X
-2- i
så är Å1,2 = -2 ± i (se exempel 3.2) och alltså är u(A) detta fall är lösningarna dämpade svängningar.
= -2. Systemet är stabilt. I □
Gränsfallet med neutralt stabila system är relativt ovanligt, men ett typiskt exempel ser ut på följande sätt.
66
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
Exempel 4.3 Matrisen w
l
u
-1
har egenvärdena Å1 = -1. Å2 = i, Å3 = -i. Alltså är u(A) = 0. Eftersom A är diagonaliserbar så är systemet neutralt stabilt. Den allmänna lösningen till (H) har formen xi(t) = acost + bsint { x2(t) = ce-t x 3 (t) = bcost - asint och är begränsad men går inte mot O då t ~ +oo om (a, b) =I (0, 0). För stora t blir lösningarna praktiskt taget rena harmoniska svängningar. D Slutligen ger vi ett instabilt exempel. Exempel 4.4 Om w
-1
så är Å1 = 1 och Å2 = -1. Alltså är u(A) allmänna lösningen till (H) är
=
l
u
1 och systemet är instabilt. Den
och dessa uttryck är obegränsade för storatom c1 =I 0.
D
4.2. FRIA SVÅNGNINGAR
67
Stabilitetskriteriet måste modifieras om vi arbetar med icke diagonaliserbara systemmatriser, men lyckligtvis endast i gränsfallet med neutral stabilitet. En icke diagonaliserbar matris har ju minst ett multipelt egenvärde Å, och ett sådant kan ge upphov till, förutom eÅt, även termer av formen teÅt, t 2 eÅt, . . . i den allmänna lösningen till (H). Men eftersom exponentialfunktionens belopp växer mot oändligheten (om Re Å > 0) eller avtar mot noll (om Re Å < 0) mycket snabbt i förhållande till variationen i potensfaktorn t•, så kan vi få en ändrad slutsats endast i gränsfallet u(A) = 0. Multipla egenvärden Å = iw på imaginära axeln kan ibland ge upphov till termer och alltså instabilitet. V1 sammanfattar:
Sats 4.2 För varje lineärt homogent system gäller u( A) < 0 u( A) > 0 u( A) = 0
==> ==> ==>
x = Ax
med konstant systemmatris
stabilitet instabilitet neutral stabilitet eller instabilitet.
Bevis: Samma resonemang som i det diagonaliserbara fallet.
D
Om u( A) = 0 och alla egenvärden på imaginära axeln är enkla, så har vi säkert neutral stabilitet, medan om det finns ett multipelt egenvärde där så är instabilitet möjlig men inte garanterad.
Exempel 4.5 Om A= [~ ~] så har vi ett dubbelt nollställe Å ningen blir
x( t)
=0
= C1
på den imaginära axeln. Den allmänna lös[
~
]
+ c2
och systemet är neutralt stabilt.
[
~
]
D
Exempel 4.6 Om A
= [~
~
så har vi återigen ett dubbelt egenvärde i Å
l =
0. Denna gång är matrisen inte diagonaliserbar. Den allmänna lösningen till (H) fås lätt, eftersom A är triangulär. Den blir
68
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
Alltså är systemet instabilt. eftersom t-tennen inte är begränsad för stora t.
□
I det stabila fallet går varje lösning mot noll för stora t, men det är också av intresse hur snabbt man kan vänta sig att detta sker. Låt oss för enkelhets skull se på det diagonaliserbara fallet. Den allmänna lösningen kan då skrivas n
x(t)
= L c1ce>..•ts1r lr=l
och en uppskattning med triangelolikheten ger n
jx(t)I::;
L lc1rlle>..•tlls1rl::;
Keu(A)t
lr=l
då t
~
0, ty
ie>..•'I = eReÅ•t
::; eu(A)t
då
t
~
0,
för varje egenvärde Å1r till matrisen A.
Sats 4.3 Om A är diagonaliserbar, så uppfyller varje lösning x(t) till (H) en olikhet
lx(t)I::;
Keu(AJt.
t
~
0.
(Konstanten K beror förotom av systemet också av lösningens begynnelsevärden.) I det stabila fallet är u(A) < 0. Vi inför då tidskonstanten för systemet genom formeln 1 Ta= - a(A)'
Ta
(a för amplitud)
Olikheten ovan kan då skrivas som
Tidskonstanten för ett stabilt system anger den tidsskala i vilken en allmän fri lösning går mot noll. Efter t ex 10 tidskonstanter har lösningens amplitud avtagit med en faktor e- 10 :::::: 10- 4 _ I stort sett gäller samma sak även för icke diagonaliserbara system, eftersom de där tillkommande faktorerna tk varierar så långsamt i förhållande till exponentialfaktorema för stora t. Man bör noga hålla i minnet att stabilitet i den ovan definierade meningen inte är en egenskap hos ett verkligt system utan hos en matematisk modell som vi gjort av ett sådant. Man bör därför vara försiktig med slutsatser från modellen till verkligheten eller omvänt. Medvetna om denna risk nämner vi i alla fall några allmänna resultat om stabilitet av modeller för praktiska system. Varje passivt RLC-nät är stabilt, eller i några undantagsfall neutralt stabilt. För aktiva nät, som innehåller
4.2. FRIA SVÅNGNINGAR
69
t ex negativa resistanser (realiserade genom komponenter som tunneldioder, transistorer eller operationsförstärkare) eller förstärkare, kan man få instabilitet. Detta torde vara bekant för var och en som försökt konstruera en förstärkare. Konservativa mekaniska system är antingen neutralt stabila eller instabila, medan icke konservativa system (tex med friktion) kan vara stabila i vår mening. Slutligen en terminologisk anmärkning. Det som vi ovan kallat stabilitet brukar mycket ofta ges namnet asymptotisk stabilitet. Skälet till att vi valt det kortare namnet är främst bekvämlighet, men också att vår inskränkta definition i många fall ger det enda praktiska stabilitetsbegreppet, medan den neutrala stabiliteten svarar mot ett speciellt gränsfall. I t ex kursen i Olineära dynamiska system används den mera invecklade terminologin, och där måste man skilja på flera olika stabilitetsbegrepp, vilka dock sammanfaller i det enkla fall som vi arbetar med i den här kursen. Vi skall nu använda teorin ovan på en mycket grov modell av ett tyvärr alltför bekant tekniskt tillämpningsexempel.
Exempel 4. 7 (Neutrondiffusion) I exempel 2.7 undersökte vi neutronkoncentrationen i en uranstav och fick fram tillståndsekvationema
Om staven består av U235 så fångas en del av neutr 0 k > h2 h >
{D
Vk
etc. Eftersom h = l/2 så ser vi att stabiliteten beror av stavens längd l. Resultatet är att systemet byter karaktär vid den kritiska längden
stabilt neutralt stabilt (»reaktor») instabilt (»bomb»). (Lyckligtvis har kärnkraftverk inte lineär dynamik av ändlig ordning, utan det finns olika stabiliserande faktorer.) D
4.3
Tvungna svängningar
Vi övergår nu till system med insignaler. Vi beskriver dem med ekvationer av typen dx dt =Ax+ f(t)
(T)
där f(t) = (/1 (t), f 2 (t), ... , fn(t)) är givna funktioner som bestäms av insignalema. För system av denna typ gäller, liksom alltid för lineära ekvationer, att den allmänna lösningen lätt kan erhållas om man känner den allmänna lösningen till motsvarande homogena problem och en enda lösning till den inhomogena. Vi har följande bekanta princip: Sats 4.4 Om xp är en lösning till det inhomogena systemet (T), så är den allmänna lösningen till detta system av formen
x(t)
= Xp(t) + XH(t),
där x8 betecknar den allmänna lösningen till det homogena systemet (H). Om speciellt A är diagonaliserbar, så är den allmänna lösningen till (T) av formen x(t)
= c1eA '5i + · · · + c;.eA"'Sn + xp(t). 1
4.4. INSIGNALSTABILITET
71
Bevis: (Se också Sparr (1995).) Summan d dt(xp
+ XH)
dxp
= dt
dxH
+ dt
= Axp
Xp
+ XH är en lösning till
(T), ty
+ f(t) + AxH = A(xp + xH) + f(t).
Att varje lösning kan skrivas på detta sätt följer av att d dt (x - xp)
=
dx dxp dt - dt =Ax+ f(t) - (Axp
+ f(t)) = A(x -
xp).
Detta visar nämligen att x - Xp är en av lösningarna XH till det homogena systemet. D
Den speciella lösningen xp brukar traditionellt kallas för en partikulärlösning. Den är dock inte mer speciell än att den kan bytas ut mot vilken som helst annan lösning till den inhomogena ekvationen. Då blir formen på den allmänna lösningen densamma, men konstanterna c1 , c2 , får andra värden. Det är nu fallet att i det långa loppet spelar det varken för stabila eller instabila system någon större roll vilken partikulärlösning som valts vid lösningen. Vi skall se att för stabila system blir alla partikulärlösningar nästan lika efter tillräckligt lång tid. För instabila system så överväger däremot de fria svängningarna för det mesta över en måttlig insignal. Endast snabbt växande insignaler kan överväga över de växande fria svängningarna. Antag först att det fria systemet x = Ax är stabilt. Då gäller för varje lösning XH till detta att XH ➔ 0 då t ➔ +oo. Lösningen till det inhomogena systemet har formen x(t)
= Xp(t) + XH(t)
Efter ett antal tidskonstanter blir alltså x( t) praktiskt taget identisk med den partikulärlösning vi valt. Alla partikulärlösningar blir ungefär lika efter tillräckligt lång tid. Vi kan uttrycka detta så att ett system med stabil systemmatris så småningom »glömmer bort» sina begynnelsevärden. När detta skett, beror lösningen x(t) bara av vilka värden insignalen f(t) haft. Hur lång tid detta tar beror i huvudsak på hur stor systemets tidskonstant är. De homogena lösningarna kallas i det stabila fallet ofta för transienter (transient = övergående).
4.4
Insignalstabilitet
För ett system med insignaler bestäms uppförandet av tillståndet x( t) för stora t inte bara begynnelsevärdena utan även av utseendet på insignalen. I detta fall är det lämpligt att tala om stabilitet i en något annan mening än för fria svängningar. Systemet sägs vara insignalstabilt om måttligt stora insignaler endast kan ge upphov till måttliga tillståndsvariationer. Vad som skall menas med måttlig signal kan
72
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
variera, beroende på ändamålet med undersökningen. För vissa ändamål använder man absolut integrerbara signaler ( L 1-signaler). för andra kvadratiskt integrerbara, eller mer fysikaliskt uttryckt, signaler med ändlig totalenergi (L2-signaler), medan det ibland är naturligast att betrakta signaler som är begränsade (för stora positiva t) som måttliga. För vår klass av systemmodeller, alltså lineära tidsinvarianta tillståndsmodeller, kan man visa att dessa typer leder till precis samma stabilitetsbegrepp, och vi gör följande definition: DEFINITIOI'\ 4.2 Ett system på tillstånds/orm sägs vara insignalstabilt om för varje begränsad insignal tillståndsvariablerna blir begränsade funktioner av tiden.
För system av formen (T), alltså lineära tidsinvarianta av ändlig ordning, så visar det sig lyckligtvis att denna nya typ av stabilitet gäller för exakt samma typ av system som har stabilitet av fria svängningar. Vi kan nämligen bevisa följande sats:
Sats 4.5 Om matrisen A är stabil och om f(t) är begränsad då t lösning till systemet dx dt =Ax+ f(t)
~ t0
så är varje
begränsad då t ~ t 0 • Om A är neutralt stabil eller instabil, så finns begränsade funktioner f(t) som ger upphov till lösningar x(t) som är obegränsade för stora t. I det instabila fallet är påståendet självklart, ty då finns obegränsade lösningar till och med när insignalen är identiskt noll. Bevisen i de stabila och neutralt stabila fallen kommer i nästa kapitel. På grund av satsen ovan kan vi alltså för system av typen (T) tala om stabilitet utan att behöva ange om det gäller stabilitet mot ändringar av begynnelsevillkor eller av insignaler. I föregående kapitel fann vi en nästan allmän metod att bestämma lösningar till systemet (T). I nästa kapitel skall vi göra en elegant omformulering och utvidgning av denna metod, som ger lösningen även dä A inte är diagonaliserbar. Här övergår vi nu i stället till en enklare metod att bestämma partikulärlösningar. Denna fungerar dock endast i vissa speciella, fast praktiskt mycket viktiga, fall.
4.5
Stationära lösningar
Om insignalen har ett exponentiellt tidsberoende, alltså av formen e•t med komplex frekvens s, sä är det naturligt att söka en partikulärlösning som beror av tiden pä samma sätt. Denna metod torde vara bekant från envariabelanalysen. Lät alltså drivtermen ha formen f (t) = e•t f, där / är en given konstant vektor, och ansätt lösningen pä formen x(t) = e•tz, där z är en sökt konstant vektor. Insättning i (T) ger ekvationen
se•tz
= Ae•tz + e•t f.
4.5. STATIONÅRA LÖSNINGAR
73
Eftersom en # 0 för alla t, så är detta ekvivalent med ett lineärt ekvationssystem för z på matrisform: (si - A)z = /. Matrisen för detta system är si - A. Som bekant uppträder vid lösningen två olika fall, beroende på om matrisen si - A är singulär eller ej: • si - A inverterbar: I detta fall finns en entydig lösning, given av formeln
z
= (si -
A)- 1 I.
• si -A singulär: I detta fall har i regel ekvationssystemet ingen lösning alls. I undantagsfall, tex alltid då/= 0, så finns det oändligt många olika lösningar.
Lägg märke till att differentialekvationen alltid har lösningar, även i det andra fallet. I detta fall har dock i regel ingen lösning ett rent exponentiellt tidsberoende. När inträffar fall 2, alltså när är si -A inte inverterbar? Ett välbekant kriterium fås ur determinantteorin. Enligt denna är en matris singulär då och endast då dess determinant är lika med noll. Men det(sl - A) = 0
~
PA(s) = 0
~
s är ett egenvärde till A.
Fall 2 inträffar alltså precis då s är ett egenvärde till A, och vi har för tredje gången träffat på egenvärdesproblemet vid våra försök att lösa differentialekvationer. De flesta tal är ju inte egenvärden till A, och lösningsansatsen lyckas alltså i allmänhet.
Sats 4.6 Om insignalen till systemet (T) har komplex frekvens s, alltså f (t) = e•t f, där s inte är en egenfrekvens till systemet, så finns en lösning till systemet som beror exponentiellt på tiden med komplex frekvens s. Denna lösning är den enda med ett sådant tidsberoende, och ges av formeln
Om insignalen är en (komplex) harmonisk svängning, alltså om s = iw är rent imaginärt, så blir den ovan funna lösningen också en harmonisk svängning, och är speciellt begränsad för alla t. Denna lösning brukar i tekniska sammanhang kallas för den stationära lösningen till systemet (T), och är mycket flitigt använd. Formellt sett finns det ingen större skillnad på den funna lösningen då s ligger på den imaginära axeln och då s ligger utanför denna. Det senare fallet kommer längre fram att leda till det viktiga begreppet överföringsfunktion. Vi inför benämningen generaliserat stationär lösning för funktionen
74
KAPITEL 4. STABILITET OCH STATIONÄRA LÖSNINGAR
I fallet 2 bar vi en insignal. vars komplexa frekvens är en egenfrekvens till systemet. Insignalen bar då samma tidsberoende e'' som en viss fri svängning, en mod, till systemet. Detta fall kallar vi för resonansfallet. ty det bar nära samband med fysikaliska resonansfenomen, speciellt då s = iw är rent imaginärt. Innan vi använder lösningsmetoden ovan på konkreta exempel är det nog lämi:r ligt att repetera rnatrisinversion.
4.6
Resolventmatrisen
I lösningsforrnlerna ovan ingår matrisen RA(s) ~ (si - A)- 1 ,
sej egenvärde till A.
Denna matris kalla.s för resolventmatrisen till A (eller till systemet). Den upi:r kommer även då man använder andra lösningsmetoder som Fouriertransformationen eller Laplacetransformationen. Vid räkning på dator är det som bekant ofta.st lämpligt att bestämma matrisinversen med någon variant av Gausselirnination eller snarare att lösa ekvationssystemet direkt med denna metod. Vid handräkning på små system (n = 2 eller 3) eller på systern med speciell struktur, samt för att studera allmänna egenskaper bos resolventen, är däremot en deterrninantforrnel för matrisinversion mycket nyttig. (Men använd den inte i Numerisk analys.) Vi påminner om denna formel. Om A = [ai,] är en n x n-matris så definieras adjunkten adj A till A på följande sätt:
On1
Stryk alltså den rad och den kolonn där a,i står. Nu gäller den viktiga formeln A adj A = (det A) I.
Om det A # 0 så följer härav att A-1
1 ad ' A = detA J .
Specialfallet n = 2 är mycket enkelt, och det är lönt mödan att lära sig tillvägagångssättet utantill. Om
4.8. RESOLVENTMATRISEN
75
Minnesregel: Byt plats på diagonalelementen och tecken på de andra. Formeln för inversen av en 2 x 2-matris blir alltså A-1 _
1
[ a22
aua22 - a21a12
-a21
-a12 ] au
Den allmänna formeln blir mycket mera invecklad redan för n = 4. Vid beräkningen av adj A så måste n 2 = 16 determinanter av ordning n - 1 = 3 tas fram. För större n blir det ännu mycket värre, så att här har formeln endast teoretiskt intresse (vilket dock inte är att förakta). Den ger väsentligen bara ett annat sätt att uttrycka utveckling av en determinant efter rad eller kolonn, och kan sägas vara en variant av Cramers regel. Vi använder nu den allmänna regeln på matrisen si -A. Eftersom det(sl -A) = PA(s) så blir resultatet
Sats 4. 7 Resolventmatrisen för en kvadratisk matris A kan beräknas med RA(s)
= (si -
A)- 1
= _!__() adj(sl PA
A).
8
Elementen i adj(sl - A) är underdetenninanter av ordning n - 1 till si - A och blir alltså polynom av grad högst n - 1 i s. Härav följer ett kvalitativt viktigt resultat:
Följdsats 4.8 Elementen i resolventmatrisen (si - A)- 1 är rationella funktioner avs, alla äkta bråk med PA(s) som nämnare. ( Äkta bråk innebär att täljarens gradtal är mindre än nämnarens.)
Exempel 4.8 Bestäm resolventen till matrisen A
= [ -~ _:
l·
Lösning: _ [ s - 1 -4 si -A 2 s +5 PA(s)
= det(sl -
adj (8 I - A)
A)
l
= s2 + 4s + 3
= [ 8 -2 +5
4 s-1
l
Alltså är R
8
_
A( ) -
s2
I [s+5 4 + 4s + 3 -2 s - I
l
D
76
4. 7
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
Generaliserat stationära lösningar och superposition
Vi kan nu återgå till differentialekvationerna. Formlerna i föregående avsnitt ger omedelbart en lösning på följande problem.
Exempel 4.9 Finn den generaliserat stationära lösningen till systemet
Lösning: Systemet är av formen ovan med
Enligt teorin är lösningen av formen
Insättning av s
= 2 i exempel 4.8 ger R (2) A
= 2_ [ 15
i 4 -2 1
l
'
Svar: Den stationära lösningen är
i}e"
x 1 (t)
=
X2(t)
= 15 e21
{
Anm Metoden fungerar även om s = 2 byts mot något annat (ev komplext) tal, D skilt från egenvärdena till A, alltså från -1 och -3. Den generaliserat stationära lösningen är ju bara en bland många lösningar till det inhomogena systemet (T). Det är naturligt att fråga sig om den har någon speciell betydelse för det praktiska systemet, utom den att den matematiska formen är enkel. Låt oss se på den allmänna lösningen. I det diagonaliserbara fallet har denna formen
Vilken blir den dominerande termen för stora t? Svaret kommer att bero på den inbördes storleken hos de reella talen ReÅ 1 , ... , ReÅn och Res. Om Res> u(A)
4.7. STATIONÅRA LÖSNINGAR OCH SUPERPOSITION
77
så blir till slut den stationära lösningen mycket större än de fria svängningarna. I motsatt fall, om Res < u(A) så dränks den i allmänhet till slut av de senare. Slutsatsen blir att om man är intresserad av vad som sker för stora tider, så är
det meningsfullt att tala om en (generaliserat) stationär lösning precis i det fall då Res > u(A), och möjligen i fallet då Res= u(A). I exempel 4.9 är den allmänna lösningen av formen x(t)
= c1e-t
[-~
l
+c2e-3t
[-! l
+ 115e2t [
!l
och alltså är x(t) :::::: Xgnat(t) för tillräckligt stora t. Drivtermer f (t) = e'' / med rent exponentiellt tidsberoende kan tyckas vara speciella, men eftersom systemet (T) är lineärt i f, så kan vi använda superposition för att finna lösningar med allmännare f.
Sats 4.9 (Superpositionsprincip Il) Betrakta systemet
dx
dt =Ax+ f.
(T)
Om x = x1 är en lösning till (T) med f = '1 och x = x2 en lösning med f = '2, så är x = c1x1 + c2x2 en lösning till {T) med f = c1'1 + c2'2-
Bevis: X = C1X1 + C2X2 ====? X = C1X1 + C2X2 A(c1x1 + c2x2) + c1'1 + c2'2 =Ax+ (c1'1 + c2'2).
=
C1(Ax1 + '1) + C2(Ax2 + '2)
= D
En analog formel gäller naturligtvis även för insignaler som är lineärkombinationer av mer än två termer. Om speciellt
f(t) där inget av talen s 1 , lärlösning av formen
x(t)
... ,
= C1e' '/1 + 1
.. . c,,e'P'/,,
s,, är en egenfrekvens till systemet, så får vi en partiku-
= c1 e' '(s 1 / 1
-A)- 1 /
1
+ ... c,,e'P'(s,,I -A)- 1 /,,
som också kan förtjäna benämningen generaliserat stationär. Det viktigaste specialfallet är insignaler med sinusformat tidsberoende. Sådana faller in under det föregående på grund av Eulers formler coswt
= ~(eiult + e-iult)
sinwt
= 21/e'"'. t
{
"-~
- e--•).
Villkoret för att den stationära lösningen skall vara dominant över transientema blir nu att u(A) < Res = 0. Detta innebär att systemet skall vara stabilt.
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
78
Exempel 4.10 Finn den stationära lösningen till systemet
Lösning: Enligt Eulers formler är f(t) = [ ~~:;:] =
elit [
Insättning i exempel 4.8 ger R (±3i)
=
tl l=
+ e- 3it
1 [ 5 ± 3i 4 ±12i - 6 -2 -1 ± 3i
A
[
1 [ 1 =f 13i -4 =f Bi 30 2 ± 4i 7 =f i
och härav följer att . 3 R•( •)
2
[
-½
1
l
1 -7- 9i = 60 [ l-3i ] .
(t)
i
2
Alltså är den stationära lösningen x
R (-3i) [ ½ A
= _!_ [ - , 60
9i 1 - 3i
l
3it
e
_½få ]-
l
= _!_ [ -7 + ~i 60 1 + 31
_!_ [ - 7 + 9i + 60 1 + 3i
l
l
l·
-3it
e
·
Eulers formler ger sedan den rella formen av den stationära lösningen. Svar: Den stationära lösningen är xi( t)
= ; 0 ( - 7 cos 3t + 9 sin 3t)
x2 ( t)
= 310 (cos 3t + 3 sin 3t).
{
Diskussionen före exemplet visar att metoden med stationär lösning alltid kan använd~ för stabila system med sinusformade insignaler. Detta resultat sammanfattar vi i följande sats, vilken utgör den teoretiska bakgrunden till många praktiska metoder.
Sats 4.10 (Växelströmslärans huvudsats) Om ett stabilt lineärt tidsinvariant system av formen (T)
dx
dt =Ax+ f.
4. 7. STATION ÅRA LÖSNINGAR OCH SUPERPOSITION
79
utsätts för sinusformade insignaler f(t)
= eiwt /
{komplex representation)
sd är varje lösning till {T) av formen
x(t)
= Xtram(t) + Xstat(t).
Tmnsienten har i det diagonaliserbam fallet formen
Xtram(t)
= C1eÅ1t5t + ... + CneÅntSn,
där c1 , . . . , Cn är konstanter som kan bestämmas ur Xstat samt begynnelsevärdena. I praktiken försvinner tmnsienten efter ett antal tidskonstanter r0 = -1/u(A), och sedan kvarstdr endast den stationäm lösningen, som är av formen
I det stationäm fallet har alla tillstdndsvariabler ett sinusformat tidsberoende med vinkelfrekvens w. Den praktiska nyttan med växelströmsläran (jw-räkning) är att det dynamiska systemet
dx
dt =Ax+ f(t)
formellt överförs i ett statiskt system (iwl - A)z
=/
vilket kan behandlas med den lineära algebrans metoder. Man kan dessutom i många fall undvika en direkt lösning av lineära ekvationssystem genom att använda från likströmsläran kända knep, såsom parallellkoppling och seriekoppling, ström- och spänningsdelning osv. Priset man får betala för att slippa derivatorna är att det blir nödvändigt att räkna med komplexa tal. Denna metod att eliminera derivator är ingalunda enbart knuten till elektrisk kretsteori, utan exakt samma tillvägagångssätt används t ex i mekanik, i optik och akustik (plana vågor) och i kvantmekaniken (reduktionen från den tidsberoende till den tidsoberoende Schrödingerekvationen). Med hjälp av superpositionsprincipen kan en mängd allmännare fall behandlas. Från Fourieranalysen är det bekant, att en tämligen godtycklig periodisk funktion kan skrivas som en summa av i regel oändligt många harmoniska svängningar, genom Fourierserieutveckling. Det är nu möjligt att utvidga superpositionsprincipen, i detta fall till summor av oändligt många termer, med hjälp av satser om termvis derivation. Vi genomför ej dessa räkningar utan nämner bara resultatet.
80
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
Om drivtermen f(t) till (T) är periodisk med period T, och om systemet är stabilt, så har det en entydigt bestämd T-periodisk lösning Uämte en mängd operiodiska lösningar). Denna lösning erhålls lättast genom Fourierserieutveckling. Om f(t) har utvecklingen +x;
f(t) =
L
c,. ( f )eilcnt
(OT = 21r)
lc=-x;
så är •OC
x(t) =
L (ikOI -
A)- 1c1c(f)ei1cnt_
k=-oo
Lägg här märke till att inget tal ikO är ett egenvärde till A, eftersom denna matris är stabil. Matrisinverserna är alltså säkert definierade. Funktionen f(t) är vektorvärd, och koefficienterna c,. = c,.(f) är alltså vektorer.
4.8
Resonans
Om drivtermen till systemet är av formen f(t) ningens amplitudvektor av formen
= e•' /,
så blir den stationära lös-
Men RA(s)
= ~() adj(s/ PA
A).
S
Om s ligger nära ett egenvärde till A, så blir PA(s) litet, och elementen i RA(s) blir (i regel) mycket stora. Detta fenomen är ett av de som i fysiken brukas förknippas med resonans.
Resonans: Om insignalens frekvens s ligger nära en egenfrekvens till systemet, så blir den stationära lösningens amplitud (i regel) stor i jämförelse med insignalens amplitud. Fysikaliskt intressant är detta i huvudsak då vi har ett stabilt system, alltså då u(A) < 0 och då insignalen är en harmonisk svängning med frekvens w, eller en lineärkombination av sådana med olika frekvenser.
4.9. INSTABILITET OCH LINEARITET
Resonans är av betydelse för relativt svagt dämpade system, där alltså u(A) är litet (och negativt). En mod är svagt dämpad om dess dämpning -u1c = - Re Å1c är liten men positiv. Resonans fås om insignalens vinkelfrekvens w ligger nära w1c = Im Å1c för en svagt dämpad mod.
81
iw
s
O'
X I
Resonans är ett av de allra viktigaste och vanligast förekommande fenomenen i fysiken och I tekniken. Detta har till följd att bestämning av egenvärden egenfrekvenser ~ resonansfrekvenser är ett mycket vanligt fenomen, från elementarpartikelfysik och kärnfysik, via byggnadsmekanik och reglerteknik, till den celesta mekaniken. Som ett exempel bland tusen kan nämnas att it ex Kalifornien och Japan är byggnadskonstruktörer i lag ålagda att bestämma egenfrekvenserna för sina projekt. Dessa frekvenser får nämligen inte ligga i frekvensbandet för jordbävningsvågor ( :'.S ca 10 Hz). I
~
4.9
*
Instabilitet och linearitet
Enligt den lineära teorin uppkommer i det instabila fallet lösningar, vars storlek växer obegränsat med tiden. I praktiken observerar vi sällan något sådant, utan lösningsamplituderna begränsas efter ett slag av någon faktor. Anledningen är i regel att de lineära modellerna, som ju ofta bygger på linearisering och alltså på att avvikelserna från ett arbetsläge antas vara små, upphör att vara en lämplig beskrivning av det verkliga systemet. För att bestämma vad som kommer att ske krävs en analys av en mer förfinad olineär modell. Detta bereder ofta stora matematiska svårigheter. Det finns inga metoder för lösning av olineära system som täcker alla fall, så som teorin för lineära tidsinvarianta system gör. Några ofta användbara hjälpmedel är geometrisk analys (fasrumsmetoder) och numerisk lösning, ofta kallad simulering. Sådana metoder behandlas i fortsättningskursen i Olineära dynamiska system. Den lineära analysen är i regel dock ett nödvändigt första steg. Enligt PoincareLiapunov gäller, att om lineariseringen av ett olineärt system kring ett jämviktsläge är stabil, så är även det ursprungliga olineära systemet stabilt nära detta jämviktsläge. Därför används lineär analys både då man vill undvika instabilitet, som vid konstruktion av reglersystem och förstärkare. och då instabilitet eftersträvas, som vid konstruktion av oscillatorer.
82
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
4.10
Reella system av andra ordningen
Det kan med fog hävdas att den väsentligaste och svåraste delen av en systemanalys (sedan väl systemet modellerats matematiskt) är att bestämma egenfrekvenserna, systemmatrisens egenvärden. Endast för andra ordningens system finns enkla formler för egenvärden och egenvektorer, eftersom den karakteristiska ekvationen då är en andragradsekvation. och alltså lätt kan lösas. Låt
Då är
= Å2 - (tr A)Å + det A. med tr A = au + a22 och det A = au a22 - a 12a2 1. Egenvärdena ges av PA(Å)
Å1.2
= ½tr A ±
J~(tr
A)2 - det A.
De är alltså reella om det A ~ (tr A)2 /4, komplexa annars. Vi kan också härleda enkla villkor för stabilitet. Sats 4.11 Om A är en 2 x 2-matris. så är A
stabil
~
tr A < 0
och
det A > 0.
Villkoren för neutral stabilitet blir lite mer invecklade: A
neutralt stabil
~
trA = 0, detA > 0 eller { tr A < 0, det A = 0 eller A=0.
I alla övriga fall är A instabil. Villkoren kan lämpligen illustreras i ett diagram över ett tr A-det 4plan.
detA
1 det A = 4(tr A)2
detA,
atabil
instabil
trA reella eg nvärden
trA instabil
instabil
4.11. REFERENSER
83
Ange som övning var i diagrammen som vi har neutral stabilitet resp multipla egenvärden. Finns det någon punkt där stabilitetstillståndet inte är entydigt bestämt? För andra ordningens system kan vi allså avgöra stabilitet (men inte alltid instabilitet) direkt ur det karakteristiska polynomet, utan att explicit bestämma egenvektorer och egenvärden. Även för högre ordningens system är detta möjligt. Vi skall studera sådana metoder längre fram i kursen, och nöjer oss därför med att nämna två enkla delresultat: a) Om A är reell och stabil (neutralt stabil) så är alla koefficienter i det karakteristiska polynomet > 0 ( ~ 0). b) En reell 3 x 3-matris A med karakteristiskt polynom PA (s)
= s 3 + as 2 + bs + c
är stabil då och endast då a > 0, b > 0, c > 0, ab - c > 0.
Läsaren uppmanas att träna sina algebraiska färdigheter genom att bevisa dessa båda resultat. Det första kan ofta användas för att visa att ett system är instabilt. däremot inte för att visa stabilitet om n > 2.
4.11
Referenser
1. Åström (1976): lärobok i (lineär) reglerteknik, där ju stabilitetsbegreppet är
grundläggande. 2. Spanne (1995b): innehåller bland annat en introduktion till stabilitet för olineära system.
84
KAPITEL 4. STABILITET OCH STATIONÅRA LÖSNINGAR
Inläsningsuppgifter: Kapitel 4 l Definiera stabilitet. neutral stabilitet och instabilitet för homogena lineära system av tillståndsekvationer. 2 Ange sambandet mellan egenvärdenas lägen och stabilitetsegenskapema för ett homogent system med diagonaliserbar systemmatris. 3 Vilka förändringar måste göras i föregående svar då systemmatrisen inte är diagonaliserbar? 4 Definiera tidskonstanten för ett system och förklara dess betydelse. 5 'C nder vilka förutsättningar på systemmatrisen kan man dra slutsatsen att om insignalen är begränsad för stora tider, så är även tillståndsvariablema begränsade?
6 Beskriv den algebraiska formeln för matrisinvers. 7 Definiera resolventmatrisen till en kvadratisk matris. Hurdana funktioner av variabeln s blir elementen i resolventen? 8 Vad menas i kompendiet med en (generaliserat) stationär lösning till ett system av tillståndsekvationer? Under vilka förutsättningar finns en entydig sådan? Härled i detta fall formeln för lösningen. 9 Formulera superpositionsprincipen för insignaler ( I I) och ange speciellt vilken betydelse den har för sinusformade insignaler. 10 Vilka är stabilitetsvillkoren för ett reellt system av andra ordningen? Tänk efter hur de kan härledas.
5 Matrisfunktioner
Ett viktigt hjälpmedel i teorin för lineära system är matrisexponentialfunktionen. Denna sammanfattar i en enda matrisfunktion alla dynamiska egenskaper hos systemet. I samband därmed studeras polynom och potensserier där variabeln i stället för ett tal är en kvadratisk matris. Ett grundläggande resultat på detta område är Cayley-Hamiltons sats.
5.1
Exponentialmatrisen
De lösningsmetoder som vi funnit för tillståndsekvationer bygger helt på kännedom om egenvärden och egenvektorer till systemmatrisen. Det finns även andra metoder, bland annat av numerisk karaktär, som inte explicit utnyttjar egenvärdena och egenvektorerna. Vi skall här formulera om diagonaliseringsmetoden på ett koncist sätt som kan användas för att ta fram sådana metoder och för andra ändamål. Systemet
{
dx dt =Ax+ f(t) x(to)
= Xo 85
KAPITEL 5. MATRISFUNKTIONER
86
= Si i ett diagonalt system
övergår efter det diagonaliserande koordinatbytet x di-1 dt
= Å1X1 + /1(t)
•
-
dx2 dt
= Å2X2 + '2(t)
.
-
-( ) dxn , _ dt = AnXn + 1n t vilket har lösningen Xt(t)
= eÅ1,(t-to)xk(to) +
l'
eÅ1,(I-T) A(r) dr.
Denna lösning kan på matrisform skrivas
I lösningsformeln uppträder nästan samma matris på två olika ställen. Vi inför en speciell beteckning för en matris med denna struktur. DEFINITION 5.1 Om D = diag(d1 , d2 , ...• dn) är en diagonalmatris sd definieras exponentialmatrisen e 0 genom
0 0
Lösningsformeln ovan kan då kortare skrivas
i(t)
= eAlx(t0) + [
eACt-T)f(r) dr
Vi skall nu återgå till de ursprungliga variablerna x. Eftersom x och f = s- 1 f så blir lösningen uttryckt i x av formen x(t)
= Si(t) = SeA(t- 41 )5- 1x(t0) + [
= Si,
i
= s-•x
SeA(t-T)5- 1 f(r)dr
där A = diag(Å 1 , Å2 , ... , Ån)- Här ser vi ännu en gång samma uttryck på två ställen, vilket inspirerar till ytterligare en definition:
5.1. EXPONENTIALMATRISEN
87
den kvadratiska matrisen A är diagonaliserbar och A = SAS- 1 med A diagonal. så definieras exponentialmatrisen eA genom
DEFINITION 5.2 Om
eA
= SeAs- 1
Vi skall strax visa att eA blir oberoende av vilken diagonaliserande matris S som vi valt. Lösningsformeln övergår nu i den ännu koncisare formen
Genom att införa exponentialmatrisen har vi fått en lösningsformel som ser exakt likadan ut som formeln i det skalära fallet. All information om systemets dynamik, om dess möjliga utveckling i tiden, är sammanfattad i matrisfunktionen e'A. Har man väl beräknat denna, så återstår ~endasb en integration för att bestämma systemets svar på en godtycklig given insignal. Om speciellt systemet rör sig fritt, alltså om /(t) = 0, så blir lösningen x(t)
= eA(t-to>x(to)
För ett fritt system beskriver alltså matrisen eA(t-tol direkt övergången från systemets tillstånd vid tiden t 0 till tillståndet vid tiden t. Några vanliga namn på eAt är systemets evolutionsoperator, övergAngsmatris eller fundamentalmatris.
Exempel 5.1 Matrisen A= diagonaliseras av
Enligt definitionen är då
och
[ -21
1 -2
l
KAPITEL 5. MATRlSFUNKTIONER
88
Lösningen till systemet
ges alltså av x(t) = e'Ax(O) = [
½e-3' + ½e-t -½e-3' + ½e-t _ !.e-3' 2
+ !.e-t 2
!.e-3' 2
+ !.e-t
l[ 2] - [ -1
2
-
Je-3' + ½e-t _ ~e-3' + !.e-t 2
l
2
Kontrollera genom insättning i såväl begynnelsevillkor som differentialekvation att detta stämmer. D Vad vi hittills gjort är bara att införa en ny beteckning. I fortsättningen skall vi • definiera eAt för alla kvadratiska matriser A. även sådana som inte är diagonaliserbara. så att lösningsformeln fortfarande gäller. • undersöka allmänna egenskaper hos e 1A. • ge alternativa metoder för beräkning av e1A, dels sådana som kan användas för numeriska räkningar. dels sådana som är lämpade för handräkning och vid teoretiska undersökningar.
5.2
Matrispotensserier
Beteckningen e1A kan upplevas som chockerade om man tänker på den ursprungliga betydelsen av exponenten; e" är produkten av n stycken faktorer e. Vi har dock redan generaliserat exponentialfunktionen ett antal gånger med lyckat resultat, senast till komplexa variabler ez. Finns det något sätt att beräkna ez som kan vara meningsfyllt även då det komplexa talet z ersätts med en matris A? I funktionsteorin har vi härlett sambandet
L Ik.1 z" = 1 + z + -21 z 00
ez =
2
k=O
1
+ - z 3 + ... 6
som gäller för alla komplexa z. Om A är en kvadratisk matris så är alla positiva potenser A" = A • A • • • A (k faktorer) definierade och vi sätter lämpligen A0 = I för alla A. Vi gör därför följande definition
5.3 Om A är en kvadratisk matris så definieras exponentia)matrisen eA genom serien 00 A ~ 1 k e = L- k! A. DEFINITION
k=O
5.2. MATRJSPOTENSSERIER
89
För att denna definition skall vara meningsfull så måste vi kontrollera. dels att serien konvergerar för varje matris A, dels att den ger det tidigare definierade eA om A är diagonaliserbar. En fullständig behandling av matrisserier grundas lämpligen på begreppen vektornonn och matrisnonn. Dessa är viktiga även i många andra sammanhang, speciellt vid analys av numeriska metoder. Det finns dock inte plats för detta här och vi nöjer oss därför här med att hantera uppkommande matrisserier formellt. För konvergensbevisen hänvisas till referenserna i kapitel 8. I övningarna finns en byggsats till några av bevisen för den som vill göra dem själv. Vi kan dock redan nu bevisa att den nya definitionen inte strider mot den gamla; om A är diagonaliserbar så ger de båda definitionerna av eA samma resultat.
Bevis: Falla) A = A = diag(d1, ... , dn)Multiplikation av diagonalmatriser sker genom multiplikation av motsvarande diagonalelement, så att
01e (även för k ~11c
L
k! D
lc=O
= D · D · · · D = diag(dt, ... , d!)
= 0). Alltså är
1. le/c 1e =~ L k! d1ag(d1,~, ... ,dn) = le=O
= diag
(f :, dt, f :! cl;, ... ,f ;, d!) le=O
le=O
_ 0 -- d"1ag (ed1 ,ed2 , ... ,ed,. ) -e.
le=O
Fall b) A = SAs- 1, A diagonal. Vi använder här en även i andra sammanhang nyttig formel:
Lemma 5.1 Om A = SBS- 1 så är A1e = SB1e s- 1 . Detta gäller för alla positiva heltal k, och om A är inverterbar även för negativa k. Det gäller nämligen att
Verifiera som övning formeln för negativa k. Detta ger
~ _!_ Ale = ~ _!_ SA1e s- 1 = S ( ~ _!_ A1e) s- 1 = L k! L k! L k! le=O
lc=O
(enligt a) = Se" s- 1 = eA.
lc=O
□
KAPITEL 5. MATRJSFUNKTIONER
90
Den allmänna teorin för matrisserier visar att serien för eA konvergerar för alla (reella eller komplexa) matriser A. För e'A får vi formeln
e'A
oc
t"
1
1
2A2 + -t3A3 + · •· = "'~ """' -A" = I + tA + -t k! 2 6 k=O
Ur denna formel kan de viktigaste egenskaperna hos exponentialfunktionen härledas:
Sats 5.1
1. e0
=I
2. e' 1 Aet2 A = e(t,-t 2 JA för alla t 1 och t2. 3 _ (e'A)-1
d
4- dt e'A
= e-tA
= Ae'A = e'AA
S. (eA) T = e(Ar) (trans'f)Onering) Den svåraste av dessa regler att bevisa är den andra. Det kan ske genom multiplikation av de båda serierna i vänsterledet. Den tredje regeln följer ur de två första genom att man sätter t 1 = t och t 2 = -t,
e'Ae-tA = e(t-t)A = eOA = I. Den fjärde regeln fås genom termvis derivation av serien,
!!_etA = !!_ ~ tlc Ale = ~ ktlc-1 A" = dt
dt~k! lc=O
L 00
=
~ lc=O
k!
tlc-1
(k - l)!A1c-1A = [j = k - 1] =
(
L
lc=l
00
ti
)
j!Ai A = e'AA.
J=O
Här kan A naturligtvis brytas ut även till vänster. Den tennvisa derivationen är tillåten, ty serien är en potensserie it med oändlig konvergensradie. Bevisa regel (5) som (lätt) övning. D Om flera olika matriser är inblandade, behöver inte de vanliga reglerna för exponentialfunktionen gälla. Anledningen till detta är att matrismultiplikation i allmänhet inte är kommutativ, alltså att AB# BA
för de flesta matriser A och 8. Detta har som följd att matriserna eA+s,
eAeB
och
e8 eA
i allmänhet är olika. Med hjälp av seriemultiplikation kan man emellertid bevisa regeln
5.2. MATRISPOTENSSERIER
91
Sats 5.2 Om matriserna A och B kommuterar, AB
= BA så är
Detta är en generalisering av räkneregeln (2) ovan. Vi har nu gett en generell definition av eA. Tyvärr kan serien endast i speciella fall användas för praktiska beräkningar, eftersom det i regel är besvärligt att bestämma potenserna Air. och sedan summera serien. För en del speciella matriser leder dock metoden snabbt till ett resultat.
Exempel 5.2 Beräkna e'A, där
Lösning:
A2 = [ och alltså är Air.
e•• =
i; ~
A'
~ ~
= 0 då k 2'.: 2.
l [~ ~ l
= [
~ ~
Definitionen ger nu
= I +tA + ~f A' + · · · = I +tA = [ ~
l
= 0'
n
+t [
~
Lägg märke till att matrisen A i detta exempel inte är diagonaliserbar.
D
Exempel 5.3 För enhetsmatrisen / gäller att
e" = e' I ty oc tlr.
e''
oc tlr.
=~ _,1r. =~-I= e'I L- k! L- k! . lr.=0
D
lr.=0
Denna enkla regel kan vi använda tillsammans med multiplikationssatsen, ty enhetsmatrisen / kommuterar med varje matris.
Exempel 5.4 Beräkna e' 8 , där B= [ 02 -32
J.
Lösning: Matrisen 8 kan splittras upp som 8 = 2/ - 3A, där A är matrisen i exempel 5.2. Eftersom / och A kommuterar, alltså IA =Al(= A), så är
e,s =e2tl-3tA =e2t1 e -3tA =e 2, 1 [ 01 -3t 1
l= [
2' e2' -3te 2t
0
e
l
-
D
KAPITEL 5. MATRJSFUNKTIONER
92
Lägg dock märke till att en sådan uppsplittring endast kan användas för kommuterande matriser.
5.3
Lösning av tillståndsekvationer
Räknereglerna för e'A gör nu det möjligt att lösa systemet
dx dt =Ax+ f(t) precis som i det skalära fallet. Inför en ny variabel genom
z(t) = e-Atx(t)
x(t) = eA'z(t).
p(A)
= Sp(B)s- 1 .
Detta följer omedelbart genom summation ur motsvarande regel för p( s) = s", och den har vi bevisat i ett tidigare avsnitt. Definitionen ovan av funktion av matris kan utvidgas från polynom till vissa potensserier (»polynom av oändlig grad») vilket vi redan sett för funktionen f (s) = et•. En sådan utvidgning är alltid möjlig om f är en hel analytisk funktion. Enligt funktionsteorin kan en sådan utvecklas i en potenserie
som konvergerar för alla s. Vi definierar då / ( A) genom formeln
=L 00
f(A)
=L 00
et1rA"
lr=O
J(") ( ) k! O A".
lr=O
Vi skall senare bevisa att denna serie konvergerar för alla kvadratiska matriser A. Om funktionerna/, g och h uppfyller h(z)
= f(z)g(z)
h(A)
= f(A)g(A).
såär
5.4. MATRlSPOLYNOM OCH CAYLEY-HAMILTONS SATS
95
Detta följer med gränsövergång ur motsvarande polynomformler. Sålunda kan vi definiera till exempel
och
Vi skall nu se att beräkningen av /(A) kan återföras på ett polynom av grad n - 1 i A, om egenvärdena till A är kända. Detta följer ur ett berömt resultat av två av algebrans och matristeorins pionjärer, Hamilton och Cayley.
Sats 5.5 ( Cayley-Hamiltons sats) Varje kvadratisk matris uppfyller sin egen karakteristiska ekvation PA(A)
= 0.
Det finns ett mycket enkelt bevis om A är diagonaliserbar, medan det allmänna fallet är något svårare.
Bevis: (I) A diagonaliserbar: A = SAs- 1, A = diag(Å1, ... , Ån) ~ PA(A) = SpA(A)s- 1 = 5 diag(pA(Åi}, ... , PA(An))s- 1 = = S diag(O, ... , O)s- 1 = 0 ty PA(Å1:)
= 0 för varje egenvärde Å1:.
D
Bevis: {Il) A allroäm Vi utgår från matrisidentiteten (si - A) adj(sl - A)
= det(sl -
A)I
= PA(s)I
Matrisen adj(sl - A) har element som är polynom i s av grad högst n - 1. Alltså är adj(sl - A) = Bo+ s81 + · · · + sn-l Bn-1 där Bo, Bi, ... , Bn-l är konstanta n x n-matriser. Denna formel ger efter multiplikation med si - A
-ABo + s(Bo -
ABi)
+ s 2(81
- AB-i)+···+ Sn-l(Bn-2 - ABn-1) + snBn-1 = = PA(s) I= ctol + Bct1I + · · · + sn-lctn-11 + snl
KAPITEL 5. MATRISFUNKTIONER
96
Detta är en identitet i s. och alla koefficienter m~te alltså vara lika. Om vi ersätter alla s i formeln med A. så gäller fortfarande likhet, eftersom på båda sidor alla s står till vänster om sina koefficienter. Detta ger
-ABo + A(Bo - ABi) + A2 (B1 - AB.i)+ ... + A"- 1 (Bn-2 - ABn-d + A"Bn-1 = = ao/ + Ao1 + · · · + A"- 1an-1l + A"I Vänsterledet blir nu en teleskopsumma med värdet noll. medan högerledet just är □
PA(A).
ur Cayley-Hamiltons sats
ser vi att A" kan uttryckas som en lineärkombination av /, A, ... , A"- 1 . Multiplikation av detta uttryck med A visar att A"~ 1 = AA" kan uttryckas med hjälp av A. A2 . ... ,A" och alltså även med /, A, ... , A"- 1 . Successivt kan sedan detsamma göras för varje potens A" av A och summation visar att varje polynom i A kan reduceras till ett polynom av grad högst n - 1. Med lite överblick inser man att denna reduktion i själva verket endast är en polynomdivision med PA och slutresultatet är resten vid denna division. Exempel 5.5 Antag att matrisen A har det karakteristiska polynomet PA(s) s 2 + 2s + 3. Skriv A4 som en lineärkombination av / och A. Lösning: Polynomdivision av p(s) = s4 med PA(s) = s 2 + 2s + 3 ger
s4
= (s 2 -
2s + l)(s 2 + 2s + 3) + 4s - 3.
Insättning av A ger då, eftersom A2 + 2A + 3/ A4
= (A2 -
2A + /)0
= 0, att
+ 4A -
3/
= 4A -
3/
vilket är det önskade uttrycket.
D
För potensserier, som i någon mening är polynom av »oändligt högt gradtal», borde en liknande reduktion gå att göra. Vi skall bevisa detta och också ange en algoritm för beräkning av restpolynomet. Lemma 5.3 Om f(s) är en hel analytisk funktion och p(s) är ett polynom av grad n, så kan f skrivas som f(s) = g(s)p(s) + q(s)
där g( s) är hel analytisk och q( s) är ett polynom av grad < n.
5.4. MATRISPOLYNOM OCH CAYLEY-HAMILTONS SATS
97
Bevis: Antag att p har nollställena Å1 , ... , Åm med multipliciteterna v1 , "2, ... , Vm respektive, alltså att p(s)
= k(s -
Å1 t1
· · · (s ~
Då har / (s) / p( s) poler av ordning
/(s) -(s) p
=L m
v1 + · · · + vm
Åm)"m,
=n
v,e i s = Å1c och kan alltså skrivas som
L c1c1(s ""
I
Å1ct
+ g(s)
lc=l l=l
där g är hel analytisk. Multiplikation med p(s) fullbordar nu beviset av lemmat. Låt nu p = PA och skriv
f(s)
□
= g(s)pA(s) + q(s).
Då är f(A)
= g(A)pA(A) + q(A) = q(A)
enligt Cayley-Hamiltons sats. Hur kan man praktiskt bestämma polynomet q(s)? En möjlighet är att använda metoden i beviset av lemmat, som enligt funktionsteorin går tillbaka till en Taylorutveckling av (s - Å1c)"" f(s)/PA(s) av ordning v1c - 1 i s = Å1c, för k = 1, ... , m. En annan möjlighet är ibland lättare vid enkel handräkning. Den bygger på att
f(s) - q(s)
= g(s)pA(s) = g(s)(s -
Åi}"1
•••
(s - Åm)"m.
Villkoret på q är att / (s) - q( s) har nollställen av ordning minst v1c i s = Å1c för alla k. Detta villkor är ekvivalent med att/ och q har samma derivator av ordning < v1c is= Å1c, k = 1, ... , m, återigen på grund av Taylors formel. Vi sammanfattar:
Polynom.metod för beräkning av / ( A). Antag att A har de olika egenvärdena Å1 , ... , Åm med multipliciteterna v1 , "2, ... , lim- Om f är hel analytisk sd finns ett polynom q av grad < n som uppfyller de sammanlagt n ekvationerna
= /(Å1c) q (Å1c) = /'(Å1c) q(Å1c) 1
q(111,-l)(Å1c)
me.d k
= 1, ... , m.
= f'"1s-l)(Å1c)
För detta q är f(A)
= q(A).
KAPITEL 5. MATRlSFUNKTIONER
98
Exempel 5.6 Beräkna med polynommetod e'A. där A
=
[ -21 -21
l
.
Lösning: Matrisen A har enkla egenvärden .Å 1 = -1 och .,\ 2 = -3. Sök ett polynom = q( A) ! Villkoren blir för / (s) = e''
q sådant att e'A
{
q(-1)
= /(-1) = e-t
q(-3) = J(-3) = e- 3t_
Alltså är
etA = q(A) = ½(3e-t - e-3t)/ + ½(e-t - e-3t)A
= ½(3e-t -
e-3t) [
l
~ ~ + ½(e-t -
e-3t) [
-i -~ l·
Jämför detta med resultatet i exempel 5.1.
D
Vi skall också se på ett besvärligare exempel med ett multipelt egenvärde.
Exempel 5. 7 Antag att A har det karakteristiska polynomet PA(s)
= (s + 2)(s -
1) 2 •
Skriv e'A som ett andragradspolynom i A! Lösning: Ansätt q(s) = a 0 + a 1 s + a 2 s 2 . Villkoren blir, med f(s) q(-2) { q(l)
q'(l)
= /(-2) = /(1) = /'(1)
ao { ao
= e"
+ 402 = e- 2t + a 1 + a 2 = et 01 + 202 = tet
201
Detta lineära ekvationssystem i ok löses till exempel med Gausselimination och resultatet blir 1 -2t 8 t 2 t ao = -e + -e - -te 9 9 3 2 -2t 2 t 1 t 01 = --e + -e + -te 9 9 3 02
1 -2t = -e 9
1 9
t
1 3
t
-e + -te .
5.5. BEVIS AV STABILITETSSATSERNA
99
Detta ger till slut
Svar: etA
= !(e- 2' + Be' -
6te')I + !(-2e- 2'
+ 2e' + 3te')A + !(e- 2'
-
e' + 3te')A2 . D
Den som genomfört räkningarna i exemplen inser att polynommetoden för etA allmänt leder till ett lineärt ekvationssystem med högerled av formen t'e>.,.t, 0 ~ l < V1t- I praktiken är detta system besvärligt att lösa och metoden i beviset av lemmat är oftast lika enkel. Lägg märke till att denna metod är nära besläktad med partialbråksuppdelning och residyberäkning. Den viktigaste behållningen av metoden är den allmänna beskrivningen av exponentialmatrisens utseende:
Sats 5.6 Varje element i exponentialmatrisen etA är en lineärkombination av termer av formen t'e>.t, där ;\ genomlöper spektrum av A och l är mindre än multipliciteten för ;\.
Beräkningsmetodema i detta avsnitt är av principiellt intresse men oftast ganska obekväma i konkreta fall. Vi kommer i samband med Laplacetransformationen att beskriva en något mer hanterlig algebraisk metod. Som en utmaning till läsaren ger vi slutligen två funktionsteoretiska varianter av denna metod. Låt RA(s) = (si - A)- 1 vara resolventmatrisen till A. Den är en rationell funktion av s med poler i egenvärdena till A. Då gäller
Cauchys integralformel för matrisfunktioner: Om C är en kurva som enkelt omsluter spektrum av A och om f är analytisk på och innanför C, så är
f(A)
= -21 . 1n
[ RA(s)f(s) ds.
le
Residyformeln för matrisfunktioner: Om f är analytisk på spektrum för A, så är Res(RA(B )f (s) ). f(A) = , egenvärde till A
5.5
Bevis av stabilitetssatserna
Vi kan nu bevisa stabilitetssatsema i förra kapitlet även i fallet då A inte är diagonaliserbar. Inför tillfälligt beteckningarna
lxl = l:5J:5n m?,X lx;I, IAI = 1:51,1:5n ~~ la.;;I
KAPITEL 5. MATRlSFUNKTIONER
100
En allmännare behandling av vektor- och matrisnonner kommer i (ett planerat) kapitel 8. En enkel räkning visar att
IAxl Eftersom
laijl
~
IAI
och
lxjl
lxl
~
~
ni Al lxl
för alla i och j är ju
n
n
IAxl = l.. - Å)z* z
= 0.
Om z # 0 så är z• z > 0 och det följer att >.. - X = 0. Alltså är varje egenvärde till Areellt. D Egenvektorerna till en symmetrisk matris A kan då också väljas reella, ty de fås som lösningar till de reella ekvationssystemen (A - >..l)z
= 0.
Ännu en viktig egenskap hos dessa egenvektorer är följande:
8.8. SYMMETRISKA MATRISER
115
Sats 6.10 Om A är en reell och symmetrisk matris så är (reella) egenvektorer som hör till olika egenvärden till A ortogonala.
Bevis: Om A är reell och symmetrisk så är A* = A. Antag att Az1 = Å1z1 och att Az2 = Å2z2. Då kan uttrycket zi7 Az2 skri~ om på två olika sätt:
och Alltså är (Å1 - Å2)zi7 z2 = zi7 Az2 - zi7 Az2 = 0 och om Å1
i- Å2 så är zi7 z2 = 0.
D
Varje vektor skild från nollvektorn kan normeras genom att den divideras med sin längd. Om A är en reell symmetrisk n x n-matris med n olika egenvärden Å1 , ... , Ån så bildar alltså motsvarande normerade egenvektorer q1 , ... , qn ett ortonormerat system, och A diagonaliseras av den ortogonala transformationsmatrisen Q = [q1 · · · qn]-
Exempel 6.2 I exempel 6.1 uppträder den symmetriska matrisen
Dess karakteristiska polynom är PA(s)
= s3 -
8s 2 + 19s - 12 = (s - l)(s - 3)(s - 4)
och egenvärdena är mycket riktigt reella, 1, 3 och 4. Motsvarande egenvektorer blir Å= 1 : Z1 = c1(l, 2, 1) 3: Z2 = C2 (1, 0, -1) 4: Z3 = C3 ( 1, -1, 1)
Ck
#- 0
och man ser genast att de är ortogonala. Genom att välja
C1 = 1(1, 2, 1)1-1 = (12 + 22 C2 = 1(1,0, -1)1-l C3 = I(1, -1, 1) 1- 1
+ 12)-1/2
1/y16 1/y'2 1/ v'3
KAPITEL 6. SYMMETRISKA MATRISER
116
så får vi en ortonormerad bas av egenvektorer
och en ortogonal transformationsmatris
Q=
1
1
1
v'6
v'2
v'3
2
0
v'6 1
1
v'6 - v'2
1
v'3 1
v'3
som diagonaliserar A. Kontrollera som övning att Q T Q 7 AQ är diagonal.
=
Q- 1 och att Q- 1 AQ
= D
För allmänna icke symmetriska matriser kan multipla egenvärden medföra att matrisen inte blir diagonaliserbar. Detta kan inte inträffa för symmetriska matriser. För dessa gäller
Sats 6.11 (Spektralsats tör symmetriska matriser) Varje reell symmetrisk matris kan diagonaliseras med en ortogonal matris Q,
Kolonnerna i Q utgör en ortonormerad bas, besttlende av egenvektorer till A.
Vi har ovan bevisat satsen i det fall då alla egenvärden är enkla. Ett bevis i fallet med multipla egenvärden kommer i appendix B. Om det finns multipla egenvärden så uppkommer hela plan av egenvektorer. Huvudaxlarna, egenvektorernas riktningar, är inte längre entydigt bestämda, men kan väljas ortogonala. Om Q diagonaliserar A, så är
A = QAQ- 1 = QAQ 7
= [ ql
'12
Å1
0
0
Å2
0
0
qn]
o...
= q[
0 0
'Il
Ån
qJ
6.6. SYMMETRJSKA MATRISER
117
Vi har alltså formlerna
Matriserna q;qJ = Proj( q;) är projektionerna på egenriktningarna till matrisen A. Verkan av den lineära avbildningen y = Ax på en vektor x kan åskådliggöras på följande sätt: 1. Dela upp x i komposanter
q1c(q[ x) längs huvudaxlarna.
2. Sträck längs axel k i skalan Å1c. 3. Summera de sträckta komposanterna.
Denna beskrivning är mycket använd i mekaniken ( med n = 2 och 3) och i andra viktiga tillämpningar. Några viktiga storheter från tekniken som beskrivs med symmetriska matriser är:
tröghetstensorn spänningstensorn töjningstensom
egenvektorer huvudtröghetsaxlar huvudspänningsaxlar huvudtöjningsaxlar
egenvärden huvudtröghetsmoment huvudspänningar huvudtöjningar
Spektralsatsen säger till exempel för spänningstensorn, att det finns tre riktningar i rummet, där skjuvspänningen på normalplanet är noll, och att dessa tre riktningar är ortogonala. Dessa resultat är tekniskt mycket viktiga och ingalunda självklara. Multipla egenvärden är mycket vanliga i dessa tillämpningar, i regel som följd av symmetriegenskaper. Om en kropp är rotationssymmetrisk kring en axel, så är denna axel en egenriktning till tröghetstensom, och hela normalplanet till axeln består också av sådana egenriktningar, alltså huvudtröghetsaxlar.
KAPITEL 6. SYMMETRISKA MATRISER
118
6. 7
Referenser
I. Spanne (1995a): innehåller mer om spektralteori för symmetriska och andra speciella matriser.
2. Sparr (1995): grundkursen i lineär algebra. 3. Strang (1988) 4. Strang (1986)
Inläsningsuppgifter: Kapitel 6 1 Definiera skalärprodukt och ortogonalitet för vektorer i R". 2 Skriv upp och härled Cauchy-SchwarzBunjakovskis olikhet för vektorer i R" . 3 Ange matrisbeskrivningarna av speglingar, ortogonalprojektioner och vridningar. 4 Definiera begreppet ortogonal matris. 5 Definiera begreppet ortogonal matris. 6 Vad kan sägas om raderna och kolonnerna i en ortogonal matris. Varför är benämningen ortogonal lätt missvisande?
7 Definiera begreppet symmetrisk matris. Vad kan sägas om egenvärdena till en sådan? Bevisa detta. 8 Antag att A är symmetrisk och att Ax = 2x. Ay = 3y. Vad kan sägas om x och y? Bevisa detta. 9 Formulera spektralsatsen för symmetriska reella matriser. Bevisa den i fallet med enbart enkla egenvärden. 10 Vad kan geometriskt sägas om egenvektorerna hörande till ett multipelt egenvärde för en reell symmetrisk matris?
7 Kvadratiska former
Den enklaste typen av skalärvärda funktioner är konstanterna, och därnäst kommer de lineära funktionerna. Naturen kan dock ofta inte beskrivas lineärt. I första hand tillgriper man då den algebraiskt enklaste typen av olineära funktioner, de kvadratiska. Många fysikaliska storheter beskrivs bäst med hjälp av kvadratiska funktioner. Hit hör energifunktioner av olika slag men också t ex töjningstensorn för en deformerbar kropp. Kvadratiska funktioner är också av grundläggande betydelse i optimeringsläran och metoderna i kapitlet används numera flitigt för filtersyntes i telekommunikationsteorin, för att ta ett exempel av tusen.
7 .1
Matrisbeskrivning av kvadratiska funktioner
Vi inleder med ett fysikaliskt exempel.
Exempel 1.1 Vi ser än en gång på systemet av massor och lineära fjädrar i exempel 6.1. k
k
m
m
m
Lägesenergin för en lineär fjäder är ½kx2 , om k är fjäderkonstanten och x fjäderns längdändring från jämviktsläget. Summering av bidragen från de fyra fjädrarna ger 119
KAPITEL 7. KVADRATISKA FORMER
120
den totala potentiella energin V
= ½2kx 12 + ½k(x2 - xi) 2 + ½k(x3 = ½k {3x 12 + 2xl + 3x3 2 - 2x1x2 -
x2) 2
+ ½2kxl
2x2x3)
F blir alltså ett homogent polynom av andra graden i variablerna x 1 , x 2 och x 3 . Beskrivningen av F kan förenklas med hjälp av matriser. Om
-1 0] 2 -1 -1 3
och så är V= ½x 7 Kx.
Även den kinetiska energin blir en andragradsfunktion, nu i hastigheterna. Den är
T
l ·2 1 ·2 1 ·2 = 2fflX1 + 2fflX2 + 2fflX3-
J matrisform kan detta skrivas
T där M = [
= 2l x·TM·X
I ; !]
= ml.
D
Funktioner som är summor av andragradstermer brukar kallas för kvadratiska former. Form är här ett gammalt namn på ett homogent polynom, dvs ett polynom där alla termer har samma totala gradtal, och benämningen har alltså inget med modem arkitektur att göra. En allmän kvadratisk form in variabler (x 1 , x 2 , ..• , Xn) kan skrivas
f(x)
= k11X1 2 + k12X1X2 + · · · + k1nX1Xn + + k21X2X1 + k22xl + · · · + k2nX2Xn +
n
=L
~jXiXj
iJ=l
där ki; inte beror av x. V ttryckt med matriser blir detta
f(x)
= x 7 Kx
7.1. MATRJSBESKRIVNING AV KVADRATISKA FUNKTIONER
121
Exempel 7 .2 Den kvadratiska formen
kan skrivas
f(x)
= [ x1
x2
] [
~
-~ ]
[ :: ]
= x T Kx.
Framställningen är ej entydig, ty x 1x 2 -termer kan fritt bytas mot x 2 x 1-termer. Sålunda duger även D
Matrisen K kan alltid väljas symmetrisk, K båda termerna kijXiXj
=
K T, eftersom i summan av de
+ k;iXjXi = (kij + k;i)XiXj
endast summan ki; + k;i är fastlagd av funktionen /. Vi förutsätter i fortsättningen att alla matriser K som beskriver kvadratiska former har valts symmetriska. Då blir matrisen K entydigt bestämd av funktionen /.
Sats 7 .1 Antag att f är en kvadratisk form som representeras av en symmetrisk matris K. Då är K entydigt bestämd av f. Bevis: Låt först x
= e; vara vektorn med 1 på plats j
och O på alla andra. Då är
f(e;) = eJ Ke; = k;;Alltså är diagonalelementen i K alltid entydigt bestämda av /. Låt nu x = ~ + e;, med i =I= j. Då ger en enkel räkning att
Härav följer att ki;
+ k;i = /(~ + e;) - /(~) - /(e;)-
Om Kär symmetrisk, så är alltså ki;
vilket fullbordar beviset.
1
= 2(1(~ + e;) -
/(~) - /(e;)), D
KAPITEL 7. KVADRATISKA FORMER
122
7 .2
Taylorapproximation
Den vanligast källan till kvadratiska former är att allmänna olineära funktioner av flera variabler approximeras med Taylorpolynom. Låt f vara en funktion av n variabler och a = (a 1 , a 2 , ... , an) en arbetspunkt. Om man vill undersöka f i punkter nära a , så är det lämpligt att använda Taylors
formel:
som gäller om J är två gånger kontinuerligt deriverbar nära a (tillhör C2 ). Inför beteckningarna
g = /'(a) = V /(a) = [ /{(a) f~(a)
/~(a) ]
och
f{~(a) f!f:i( a)
f{~(a) f~~(a)
f:1(a) f:2(a)
f:n(a)
f:'1 (a) /~'1 (a)
K = f"(a) =
där
t = J
8f
och
OXj
f" -
å2 f
jlc -
OXjOX1c
Vektorn g är gradienten för f i punkten a. Vi definierar gradienten som en radvektor, vilket är det mest naturliga i de flesta sammanhang. Andraderivatmatrisen K kallas även för Hessematrisen (Hessianen) för f i a. Enligt ett viktigt resultat ur flervariabelanalysen är Hessematrisen alltid symmetrisk, om f tillhör C2 . Taylors formel kan med ovanstående beteckningar kompakt skrivas
Vi har skrivit funktionen 1. konstanten
f som en summa av
f (a)
2. den lineära fonnen gx 3. den kvadratiska fonnen ½x 7 Kx och resttennen termerna.
o(lxl 2 ), som nära x
= 0 är av mindre storleksordning än de övriga
7.2. TAYLORAPPROXIMATION
123
Nära en punkt a där gradienten g inte är noll. bestäms variationen av funktionen / i huvudsak av förstagradstermen f(a
+ x) - /(a)-;;::; gx för små x.
Om däremot a är en stationär punkt, alltså en punkt där f'(a) = g = 0. så är det nödvändigt att ta med ytterligare en term i utvecklingen och vi får då
f(a + x) - /(a)-;;::; ½x r Kx för små x. Uppförandet av en funktion nära en stationär punkt bestäms oftast av den kvadratiska formen x T Kx. Detta är som bekant grundläggande bland annat för optimeringsläran. Vi återkommer senare till denna men diskuterar först en situation i mekaniken där Taylorutvecklingen brukar användas.
Exempel 7.3 (Jämviktslägen för konservativa mekaniska system) Låt (q1 , 'h, ... , Qn) vara (generaliserade) lägeskoordinater för ett mekaniskt system. Om systemet är konservativt så finns en potentialfunktion V = V (q1 , 'h, ... , Qn) som anger lägesenergin i ett visst läge. De 'inre' (generaliserade) krafter som råder i systemet i konfigurationen q ges av -f(q)
= [ -v;(q)
-v;(q)
-V~(q) ]
= -VV(q)
och de 'yttre' krafter som krävs för att hålla det i jämvikt är f(q) Potentialförändringen vid en förflyttning från q = a till q = a + x är
V(a+x)-V(a) om som ovan f = [V/(a)] och K = f(a
VV(q).
= fx+ ½xrKx+o(lxl 2 )
[i,:1(a)].
Kraftförändringen blir
+ x) - f(a) = Kx + o(lxl)
Om q = a är ett jämviktsläge så är f = f (a) = 0 och alltså potentialen V stationär i a. Nära en jämviktspunkt a gäller alltså V(a + x) - V(a)
f(a + x) - f(a)
= ½x r Kx + o(lxl 2 ) = Kx + o(lxl)
Om V har ett strikt lokalt minimum i q = a så är jämviktsläget (enligt Dirichlet) stabilt. (Med våra beteckningar från kapitel 4 skulle vi säga att det är neutral stabilt, och konservativa system kan aldrig vara (asymptotiskt) stabila). I de flesta fall kan stabiliteten av jämviktspunkten avgöras genom en undersökning av den kvadratiska formen x T Kx.
KAPITEL 7. KVADRATISKA FORMER
124
För små avvikelser x undersöker man i allmänhet motsvarande lineariserade system med V ( a + x) - V ( a) = ½x 7 Kx f ( a + x) - f (a) = Kx vars rörelse i många fall är en god approximation till det ursprungliga systemets rörelse. Systemet i exempel 7.1 har förmodligen uppkommit genom en linearisering av denna typ. Matrisen K kallas för systemets styvhetsmatris, eftersom den bestämmer hur stora krafter f som behövs för att åstadkomma en given deformation x från jämviktsläget. Observera alltså att styvhetsmatrisen för ett konservativt system alltid är symmetrisk. ett viktigt resultat som omväxlande brukar tillskrivas Maxwell, Betti och Clebsch. Formeln för den potentiella energin skrivs ofta
Denna formel gäller approximativt nära ett jämviktsläge och generellt då sambandet mellan deformationer och krafter är lineärt, f = Kx. I fortsättningen av kapitlet skall vi utveckla metoder för att undersöka kvadratiska former som kan användas t ex för att avgöra om en belastad struktur, till exempel en byggnad, är stabil och för att bestämma dess knäckningslast. Samma □ metoder skall vi också använda för att lösa dynamiska problem.
7 .3
Diagonalisering av kvadratiska former
En allmän kvadratisk form är svår att analysera på grund av att den innehåller så många termer med flera olika variabler inblandade. Genom väl valda koordinatbyten kan kvadratiska former, liksom tidigare lineära avbildningar, förenklas väsentligt. Låt Vid variabelbytet övergår uttrycket för / i
där Detta är en transformationsformel som skiljer sig från den vi tidigare funnit för lineära system. Jämför
7.3. DIAGONALISERING AV KVADRATISKA FORMER
kolonnvektor x= SI I=
lineär avbildning A = sAs- 1
kvadratisk form K = (S 7 )- 1 Rs- 1
A = s- 1 AS
R = S 7 KS
s- 1 x
125
Matrisen för en lineär avbildning y = Ax transformeras vid ett koordinatbyte x SI genom en likformighetstransformation
=
A = s- 1 As medan matrisen för en kvadratisk form x 7 Kx transformeras genom en kongruenstransformation Allmänbildande anmärkning: 'Iransformationsformlerna ovan innebär att storheterna ovan är tensorer av tre olika typer. Om man inskränker sig till ortonormerade koordinatsystem så förekommer bara ortonormerade koordinatbyten och transformationsmatrisen S är ortogonal, s- 1 = 5 7 . Då bortfaller skillnaden mellan transformationslagarna för lineära avbildningar och för kvadratiska former. Man säger att de är cartesiska tensorer av samma typ.
Om
R är diagonal, S 7 KS =
R = D = diag(d1, d2, ... , dn)
så säger vi att S diagonaliserar den kvadratiska formen x 7 Kx. I de nya vari-
ablerna får vi ett enkelt uttryck
f(x)
= 1 7 D1 = d1Xi + d2x~ + · · · + dnx!.
Vi skall se att man har mycket större valfrihet beträffande transformationsmatrisen S vid diagonalisering av kvadratiska former än vid diagonalisering av lineära avbildningar. I det senare fallet måste man ju välja kolonnerna i S som egenvektorer till A. För att diagonalisera en kvadratisk form finns det däremot flera praktiska metoder att välja mellan. Den första metoden som vi tar fram återför diagonaliseringen på egenvärdesproblemet för den symmetriska matrisen K. Spektralsatsen för symmetriska matriser kan tillämpas vilket ger
Sats 7.2 Varje kvadratisk form f(x) merat koordinatbyte x blir då
=
QR ,
Q7 KQ
Q7
=
= x 7 Kx Q- 1 .
kan diagonaliseras med ett ortonorDen transformerade diagonalmatrisen
= A = diag(Å1, Å2, ... , Ån)
och den kvadratiska fonnen blir uttryckt i de nya koordinaterna f(x)
= R7 AR= Å1.ri + Å2.r~ + · · · + Ån.r!.
KAPITEL 7. KVADRATISKA FORMER
126
Kolonnerna i Q är egenvektorer till K och ,\ 1 , ... , ,\n motsvarande egenvärden.
Exempel 7 .4 Den kvadratiska formen f(x)
= 3x1 2 + 2xl + 3xl -
2x1x2 - 2x2x3
från exempel i.I överförs genom det ortonormerade koordinatbytet
x
= QR. Q=
1
1
1
v'6
J2
v'3 1 v'3
2
v'6
0 -
1
1
1
v'6
J2
v'3
på diagonalformen D
7 .4
Klassifikation av kvadratiska former Om funktionen f är stationär i a så är f(a
+ x)
- f(a)
1
= 2 x T Kx + o(lxl 2 )
där K är andraderivatmatrisen till f i a. Man kan visa att åtminstone om K inte är singulär så uppför sig f kvalitativt som den kvadratiska formen nära a, och resttermen kan för de flesta ändamål försummas. Det gäller allmänt att
• f har ett minimum i a om x T Kx > 0 för alla x :/= 0. • f har ett maximum i a om x T Kx < 0 för alla x :/= 0. • f har en sadelpunkt i a om x T Kx antar både positiva och negativa värden. Detta är enkelt att bevisa och vi hänvisar till analyskurserna. I det första fallet kallas K positivt definit, i det andra negativt definit och i det tredje indefinit. Det finns ytterligare två logiska möjligheter, nämligen • x T Kx ~ 0 för alla x, men är lika med O för något x :/= 0. • x T Kx ~ 0 för alla x, men är lika med O för något x :/= 0.
I dessa fall sägs K vara positivt respektive negativt semidefinit. Om K är semidefinit så räcker inte Taylorutvecklingen av andra ordningen till för att avgöra om f har ett extremvärde i a.
7.4. KLASSIFIKATION AV KVADRATISKA FORMER
127
Exempel 7.5 Funktionerna
och
= X1 2 -
g(X1, X2)
X~
har samma andraderivatmatris i origo, eftersom
och g
"( Xl,X 2 .) -_ [ 9'{1 g!i1 9r2 g'~
l -[ -
o 02 -12xi2
l -[ o l -
02 0
då x = (O, 0). Denna andraderivatmatris är positivt semidefinit . Funktionen uppenbarligen ett minimum i (0, 0), medan g har en sadelpunkt.
f har D
Det är alltså av praktiskt intresse att kunna avgöra om en kvadratisk form är positivt definit eller ej. Om formen är diagonaliserad, så är det mycket enkelt att avgöra vilken av typerna ovan som formen tillhör. Antag att vi har en diagonalisering f(x)
= xTKx = iTf(i = iToi = d1x~ + d2x~ + · · · + dnx!.
Eftersom koordinaterna .i- 1 , ... , Xn varierar fritt och oberoende av varandra, så ser vi att • K är positivt definit om alla di är (strängt) positiva. • K är negativt definit om alla di är (strängt) negativa.
• K är indefinit om något di är större än 0 och något är mindre.
• Kär positivt semidefinit om alla di är större eller lika med 0 och om något är lika med 0. • K är negativt semidefinit om alla
c4 är mindre eller lika med O och om något
är lika med 0. Speciellt kan man använda den diagonalisering som fås med hjälp av egenvektorerna till K, varvid d1 , ... , dn blir egenvärdena Å1 , ... , Ån- En bestämning av alla egenvärden till K kräver emellertid ett stort (och därmed dyrbart) numeriskt arbete, om antalet variabler är stort. Det är därför viktigt att det finns enklare metoder att bestämma tecknen på 0. Detta beror på att då måste f 1(x 2, ... , Xn) vara en positivt definit form i (x 2, ... , Xn) och h(x 3, ... , Xn) en positivt definit form i (x 3, ... , Xn) ... Men i en positivt definit matris måste alla diagonalelement vara > 0. Bevisa detta som (lätt) övning. Om något element d1 under räkningarnas gång visar sig vara lika med noll så kan alltså förfarandet avbrytas med svaret att K ej är definit. Vi illustrerar kvadratkompletteringsförfarandet på vår första kvadratiska form.
Exempel 7.6 Kvadratkomplettera
f(x)
= x T Kx = 3xi + 2x~ + 3x~ -
2x1x2 - 2x2X3.
Lösning: De successiva stegen blir f(x)
= 3xi -
2x 1x2
+ 2x~ -
2x2X3
1 = 3(x12 + 2(-3)X1X2) + 2x22 -
+ 3x~ 2X2X3
+ 3x32
= 3(x1 - 31 x2) 2 - 3(- 31 x2) 2 + 2x 22 - 2x2X3 + 3x32 5 3 2 5 3 2 2 = 3(x1 - 31 x2) 2 + 3 (x2 - 5x3) - 3 (- 5x3) + 3x3 = 3(x 1 - 31 x2) 2 +
5 3 2 12 2 3 (x2 - 5x3) + 5 x 3.
Vi har alltså faktoriserat K
med
R= [
i -n 1 -3 1 0
= R7 DR och
D=[~
0 5
3 0
J]
KAPITEL 7. KVADRATISKA FORMER
130
Diagonalelementen i D är alla positiva och K är alltså positivt definit. Egenvärdena till Kär ju alla positiva. Men lägg märke till att diagonalelementen i D inte är lika med egenvärdena till K. Kvadratkomplettering är inte en billig metod att beräkna egenvärden. D Räkningarna ovan ser kanske tungrodda ut, men de underlättas väsentligt om man utför dem i matrisform. I själva verket är de (vilket först konstaterades av Jacobi) identiska med de räkningar man utför vid Gausselimination, eller alternativt vid LR-faktorisering av K. Detta inses lättast genom en direkt jämförelse av beräkningsstegen: a) Dividera med diagonalelementet d 1
= k11 .
b) Eliminera så att endast nollor står under diagonalelementet i kolonn 1. c) Upprepa proceduren med kolonn 2, 3.... , n. Resultatet blir då matrisen R och de utdividerade faktorerna är diagonalelementen i D.
Exempel 7. 7 Räkningarna i exempel 7.6 blir i matrisform 1
[
-½
b) (2) + (1) 3 -1 0 ] a) d, = 3 [ 1 -1 2 -1 ---t -1 2 -10] --t 0 -1 3 0 -1 3
= 5/3 --t
a) d2
u
1
-3 1 -1
-i]
b) (3)
+ (2)
--t
[~
1
-3 1
0
-q 12
5
[~ -:] ~ -3 5
3 -1
= 12/5 --t
a) d3
u -n 1
-3 1 0
Vi kan sammanfatta resultaten på följande sätt: Sats 7 .3 Följande villkor är ekvivalenta för reella symmetriska matriser K:
1. K är positivt definit. 2. x T Kx > 0 för alla x i- 0. 3. Alla egenvärden till K är positiva. ,4. Alla pivåelement vid Gau.sselimination i K (utan rad- eller kolonnbyten) är
positiva.
7.5. KVADRATKOMPLETTERING
131
5. K = RT DR. där R är högertriangulär med ettor på huwddiagonalen, och D är diagonal med positiva diagonalelement. Det är värt att notera att faktoriseringsmetoden också är den som används för att numeriskt beräkna determinanter. Det gäller nämligen att det K = det(Rr DR)= det Rr det Odet R = (det R) 2 det D och det R = 1 medan det D
= d1d2 ... lin-
Alltså är
Motsvarande metod, Gausselimination eller ekvivalent LR-faktorisering, används även för icke symmetriska matriser. Någon eftertanke visar att metoden även ger vissa underdeterminanter till K. Följande formler gäller ku lk11 k22 k12 k21
ku k12 k13 k21 k22 k23 k31 k32 k33
= I= =
d1 d i d2
d1d2d3
och så vidare. Kontrollera i exemplet. Detta leder till följande kriterium: Den symmetriska matrisen K är positivt definit då och endast då alla determinanterna k11 k12 ku k21 k22 k23 , ... , det K ku' k31 k32 k33
är positiva. Detta kriterium anförs ofta i tillämpningslitteratur men räknemässigt är det naturligtvis klokare att faktorisera (och på så sätt samtidigt beräkna determinanterna). Faktoriseringen K = RT DR modifieras ofta något för positivt definita matriser K. I detta fall är alla d1r > 0, och vi kan dra 'kvadratroten' ur diagonalmatrisen D:
p; D= D~
2
om
D1= 2
0 0
0 ✓ p, och alltså är p = p'. □ Varje kvadratisk form kan diagonaliseras med ortogonala matriser, varvid diagonalelementen blir matrisens egenvärden. Detta visar att om K är symmetrisk och ST KS är diagonal så är • antalet positiva diagonalelement
= antalet positiva egenvärden
• antalet negativa diagonalelement
= antalet
negativa egenvärden
Antag att K och R är kongruenta, R = STKS. Om QTRQ = Dsåär (SQ)TK(SQ) = D så att K och R är båda kongruenta med samma diagonalmatris. Detta visar
KAPITEL 7. KVADRATISKA FORMER
134
Följdsats 7 .5 Antalet positiva egenvärden och antalet negativa egenvärden till en matris bevaras vid kongruenstransformationer av denna. Det enklaste sättet att reducera en kvadratisk form på diagonalform är med kvadratkomplettering. Diagonalelementen i den diagonaliserade formen blir då precis pivåelementen vid Gausselimination (utan rad- eller kolonnbyten), och deras tecken ger alltså en enkel metod att avgöra antalet positiva och negativa egenvärden till en symmetrisk matris K. Allmännare kan vi se på K - u I, vars egenvärden är egenvärdena till K minskade med u. Vi får alltså
Sats 7.6 (Spektrum.klyvning) Antag att K är en symmetrisk matris och u ett reellt tal. Då är antalet egenvärden till K som är större (mindre) än u lika med antalet positiva (negativa) pivåelement vid Gausselimination i K - ul. En förutsättning är att inget pivåelement blir 0, men detta kan bara inträffa för vissa (högst n(n - 1)/2 olika förutom egenvärdena till K) värden på u.
Övning: Tänk efter vilka värden på u det rör sig om. Med spektrumklyvning kan man enkelt avgöra hur många egenvärden som ligger i ett visst intervall. Detta kan användas som utgångspunkt för numeriska metoder för egenvärdesbestämning. En annan tillämpning är att undersöka om några egenfrekvenser ligger i ,.farliga,. frekvensband, där stora yttre störningar kan förväntas. Exempel 7 .8 Avgör hur många egenvärden till matrisen
som ligger i intervallet 2.5 < Å < 3.5. Lösning: Gausselimination i K - 2.51 ger K - 2.5/
=
0.5 -1 0 ] [ -1 -0.5 -1 0 -1 0.5
➔
[ 1 -2 0 ] 0 -2.5 -1 0 -1 0.5
➔
[ 01 -2 0 ] . 1 0.4 0 0 0.9
Pivåelementen blir alltså d 1 = 0.5, d 2 = -2.5 och d3 = 0.9, alltså två+ och ett -. Detta medför att två av egenvärdena är större än 2.5 och ett är mindre. För att se hur många egenvärden som är mindre än 3.5 Gausseliminerar vi i K - 3.5/
K - 3.5/
=
-0.5 -1 [ ~1 -1.5 -1
0 ] -1 -0.5
➔
[1
2 0 ] 0 0.5 -1 0 -1 -0, 5
➔
[1
0 ] 0 21 -2 0 0 -2.5
7.7. ODÅMPADE LINEÅRA SVÅNGNINGAR
135
Pivåelementen blir nu d1 = -0.5, d 2 = 0.5 och d3 -2.5. Antalet egenvärden över 3.5 är ett och antalet under är två. Mellan 2.5 och 3.5 ligger alltså 2 - 1 = ett egenvärde, under 2.5 ett egenvärde och över 2.5 ett egenvärde. Jämför med de exakta egenvärden vi beräknat tidigare. D
7. 7
Odämpade lineära svängningar
Nära ett jämviktsläge för ett mekaniskt system. till exempel en molekyl eller en gäller approximativt följande uttryck för den kinetiska energin T och den potentiella energin V T -- !2 x· TM·X { V= ½x 7 Kx skyskrapa,
där x = (x 1 , x 2 , ••. , Xn) betecknar avvikelserna från jämviktsläget. Massmatrisen (»tröghetsmatrisen») M och styvhetsmatrisen K är symmetriska matriser, och M är positivt definit . Enligt Lagrange blir rörelseekvationerna för systemet av formen
Mic + Kx
(L)
= g(t)
där g = (g 1 , 92, ... , 9n) representerar t ex friktionskrafter och yttre krafter på systemet. Vi hänvisar till mekanikkurserna för härledningen.
Exempel 7.9 Vi ser än en gång på systemet från exempel 6.1 och 7.1. För fria svängningar blir rörelseekvationerna mi1 { m~2 mx3
= /1 = -2kx1 = h = -k(x2 = h = -k(x3 -
k(x1 - x2)
xi) - k(x2 - x3) x2) - 2kx3.
I detta fall kan de lätt härledas genom friläggning av de tre massorna. I matrisfonn kan systemet skrivas Mic+ Kx = 0 med K
=
3k -k O ] [ -k 2k -k 0 -k 3k
D
Vi ger nu ett enklare exempel, som man med lätthet kan göra en verklig modell av och studera.
136
KAPITEL 7. KVADRATISKA FORMER
Exempel 7.10 (Kopplade pendlar) Två lika pendlar kopplas till varandra. med en kraft som beror av deras inbördes avstånd och som är noll då de båda hänger rakt ner. Hur kan de röra sig? Vi förutsätter att utslagen är små för båda pendlarna. De okopplade pendlarna har rörelseekvationerna {
ml~1 = -mg s~n 01 ml02
::::::
k
mg
-mg01
mg
= -mgsm02:::::: -mg02
Fjäderns sträckning från jämviktsläget blir :::::: 1(01 - 02 ) och fjäderkraften på den vänstra massan blir:::::: -k · sträckningen :::::: -kl(01 - 02 ), där kär Hookes konstant för fjädern. De approximativa likheterna ovan innebär alla (informella) lineariseringar av mer exakta olineära samband. De lineariserade rörelseekvationerna för det kopplade systemet blir ml~1 = -mg01 - kl(01 - 82) { ml02 = -mg02 - kl(02 - Bi) Om vi väljer vinklarna som lägeskoordinater, x 1 vationerna skrivas Mx+ Kx = 0 med
M=[ml 0
O
ml
l
'
= 01 och x 2 = 02 så kan rörelseek-
K _ [ mg + kl -kl ] -kl mg+kl .
D
Ekvationer av formen ovan kan också uppkomma som modeller för odämpade elektriska kretsar (vilka dock är ganska sällsynta fenomen):
Exempel 7.11 (Kopplade resonanskretsar) De induktivt kopplade resonanskretsarna i figuren kan beskrivas av ekvationerna L [ M där x
M ] .. [ 1/C O ] L x+ 0 1/C x
=O
= (i 1 , i 2 }. Visa detta som övning i kretsteori.
~b
c~~c D
Ekvationen (L) är inte på tillståndsform, och variablerna x = (x1 , x 2 , •.• , Xn) är ej tillståndsvariabler. Man kan överföra (L) på tillståndsform genom att välja läges-
7.7. ODÄMPADE LINEÅRA SVÅNGNINGAR
137
och hastighetsvariabler (x 1 • x 2 , ... , Xn) och (.:i: 1 , .:i: 2 , ...• .:i:n) som tillståndsvariabler. Då blir systemet (L) ekvivalent med
{!; =x
d; = -M- Kx + M-
d"
1
1g.
Detta är ett system av ordning 2n, där n är antalet frihetsgrader ( = oberoende lägeskoordinater). Det kan lösas med metoderna i föregående kapitel, men vid reduktionen förlorar man flera goda egenskaper, framför allt symmetrin hos systemmatrisen. Även om M och K båda är symmetriska, och därmed M- 1 • så är produkten M- 1 K i regel inte symmetrisk. Det visar sig därför vara lämpligare att angripa systemet (L) direkt. Härvid kan vi använda precis samma tillvägagångssätt som för lösning av tillståndsekvationer. 7. 7 .1 Separation av variabler (för den homogena ekvationen) Ansätt en synkron lösning x(t)
= .p(t)z
Insättning i (L) ger ipMz
+ .pKz = 0
och som tidigare får vi separerade ekvationer {
ip + Åip (K - ÅM)z
=0 =0
Vi leds nu till ett generaliserat egenvärdesproblem Kz
= ÅMz,
z-::/= 0
Icke triviala lösningar z # 0 fås om den karakteristiska ekvationen det(ÅM - K)
=0
är uppfylld. Denna ekvation är en polynomekvation av grad ni Å, ty utveckling av determinanten ger det(ÅM - K)
= (det M)Å" + · · · + (-1)" det K
och det M > 0 eftersom M är positivt definit . Lösningar till den karakteristiska ekvationen det(ÅM - K) = 0 kallar vi för egenvärden till paret (M, K) och motsvarande lösningar z #; 0 till Kz = ÅMz för egenvektorer eller modvektorer.
KAPITEL 7. KVADRATISKA FORMER
138
Om z är en modvektor så är z• Kz
= Az* Mz och alltså z*Kz
A=-z•Mz
Allmänt kallar man uttrycket R(K. M, z)
z*Kz
= -M z• z
för en Rayleighkvot. Denna spelar en viktig roll vid numeriska metoder för beräkning av egenvärden till symmetriska matriser. Eftersom M är positivt definit så är nämnaren z• Mz > 0 för alla z i=- 0, och eftersom M och K är reella och symmetriska så är både nämnare och täljare reella tal för alla (eventuellt komplexa) vektorer z. Detta visar att alla egenvärden är reella. Differentialekvationen ip + A'{) = 0 för amplituden I{) har olika typer av lösningar beroende på tecknet på Å. Vi får dela upp i tre olika fall: 1. Å = ÅA:
= . ..:f > 0 ger l{)A:(t) = ak cos(wkt + 61r)
2. Å = ÅA: = 0 ger ..PA:(t)
.
= ak + b1rt.
Positiva egenvärden Å1r ger alltså upphov till begränsade lösningar x(t) = 1{)1r(t)s1r medan de flesta lösningar som hör till ett egenvärde Å1r ~ 0 är obegränsade. Vi övergår nu till nästa lösningsmetod. 7. 7 .2
Diagonalisering Vi skall söka ett variabelbyte x = Si som förenklar rörelseekvati~nerna. Ur sambandet x = Si för lägena fås genom derivation ett samband x = Si för hastigheterna. Detta ger uttrycken för energifunktionerna i de nya variablerna:
Uttryckt i
x kan
rörelseekvationerna skrivas
där fA = STMS, f'< = STKS och t = ST g. Variablerna (x 1 , x2 , ... , Xn) blir frikopplade om både fA och f'< är diagonala. Vi ställs alltså inför problemet att samtidigt diagonalisera två kvadratiska former med matriser M och K. Detta problem kan lyckligtvis lösas genom att vi kombinerar tidigare funna metoder.
7.7. ODÅMPADE LINEÅRA SVÅNGNINGAR
139
Sats 7.7 Antag att M är positivt definit och att K är s-ymmetrisk. Då finns en inverterbar matris S sådan att
Kolonnerna (Si, Si, ... , Sn) i S är lösningar till det generaliserade egenvärdesproblemet för paret (M, K) och (Å 1 , Å2 , ... , Ån) är motsvarande egenvärden.
Bevis: Vi skall ta fram S i två steg. Först diagonaliserar vi M och sedan diagonaliserar vi även K utan att förstöra diagonalformen av M. Steg 1: Faktorisera M som M = er C. För att finna C kan man använda antingen Cholesky, vilket ger ett högertriangulärt Celler också spektralsatsen. vilken ger faktoriseingen M = QT OQ = (0 112 Q) T 0 112 Q, där Q är ortogonal och O är diagonal. Sätt nu S' = c- 1 . Då är
och K övergår i den symmetriska matrisen
k = S'T KS'. Steg 2: Enligt spektralsatsen finns ett ortogonalt Q' sådant att Q' T K Q' diagonal. Sätt S = S'Q'. Då är
5T MS
= A är
= Q'T S'T MS'Q' = Q'T IQ' = Q'T Q' = I
och
5T KS = Q'T S'T KS'Q'
= Q'T KQ' =
A
och det S vi har konstruerat uppfyller villkoren. Ekvationerna ST MS = I och 5T KS = A kan på kolonnform skrivas
vilket medför att ST(K - Å;M)s; följer härav att
= 0 för j = 1, ... , n.
Eftersom
sr är inverterbar
så
(K - Å;M)s;
= 0,
j
= 1, ... , n.
D
140
KAPITEL 7. KVADRATISKA FORMER
Matrisen Skallas för modmatrisen för systemet. eftersom dess kolonner är modvektorer. De mod vektorer vi tagit fram i satsen är speciella. De uppfyller ortogonalitetsrelationerna för modvektorer:
{
sI Ms1c
= 611c
sI Ks1c
= Å1cöj1c
ty [S 7 MS]j1c = s,7 Ms1c. Liksom för det vanliga egenvärdesproblemet använder man i praktiken numeriska metoder för att bestämma egenvärden och modvektorer. Den vanligaste metoden är precis den vi använt i beviset av sats 7. 7. I våra enkla övningsexempel kan det vara lättare att skriva upp och lösa den karakteristiska ekvationen, men sådana fall stöter man inte så ofta på i praktiska sammanhang. Om modmatrisen S används som transfonnationsmatris, x = Si, så övergår systemet (L) i de frikopplade ekvationerna
i~+ A1cx1c = g1c(t),
1 ::; k ::; n
Dessa ekvationer kan för givna 91c lösas med olika metoder, som behandlas i denna och tidigare kurser. I det homogena fallet, då g(t) = 0, kan vi skriva upp den allmänna lösningen uttryckt i modvektorer. Den är av formen n
x(t)
=L
cp1c(t)s1c
lc=l
där 'P1c är godtyckliga lösningar till 0. Detta är ekvivalent med att matrisen Kär positivt definit. Vi har alltså visat följande praktiskt viktiga resultat.
Sats 7 .8 Det lineära systemet
Mx+ Kx = 0
(M positivt definit)
är neutralt stabilt då och endast då K är positivt definit. Det har då n stycken egenfrekvenser (w1 , w2 , ... , Wn) som är lösningar till den karakteristiska ekvationen
det(K - w 2 M)
= 0.
Varje lösning till systemet är av fonnen n
x(t)
= L a1c cos(w1ct + 61c)s1c lc=l
1.1. ODÅMPADE LINEÅRA SVÅNGNINGAR
141
där modvektorerna St uppfyller
Ur Dirichlets stabilitetssats (se kursen i analytisk mekanik eller i olineära dynamiska system) följer att om lineariseringen av ett konservativt system kring en jämviktspunkt är (neutralt) stabil, så är jämviktspunkten (neutralt) stabil även för det ursprungliga olineära systemet. Nära jämviktspunkten är också det lineära systemets rörelse en god approximation av det olineäras. De speciella lösningarna
kallas för principalsvängningar eller normalsvängningsmoder. I en sådan utför varje lägeskoordinat en harmonisk svängningsrörelse, och alla koordinater svänger med samma vinkelfrekvens, en egenfrekvens till systemet. Beloppet av modvektorns komponenter anger det inbördes förhållandet mellan de olika koordinaternas amplituder, och tecknen om de svänger i fas eller motfas. Den allmänna lösningsformeln visar att varje fri svängning till ett lineärt system är en superposition av moder. Ett system med rörelseekvationen
kallas för en harmonisk oscillator och vi har ovan visat ett för fysiken mycket viktigt resultat: Varje lineärt (neutralt) stabilt konservativt system med n frihetsgrad.er är ekvivalent med n oberoende harmoniska oscillatorer, vars frekvenser är systemets egenfrekvenser.
Detta resultat kan utvidgas till många svängande system med oändligt många frihetsgrader, t ex vibrationer av deformerbara fasta kroppar, vätskor eller gaser, eller elektromagnetiska fält, och är också grundläggande för kvantteorin för fält. Vi tillämpar nu teorin på ett exempel.
Exempel 7 .12 Vi skall bestämma den allmänna rörelsen för systemet i exempel 7.9. Systemmatriserna
KAPITEL 7. KVADRATISKA FORMER
142
ger den karakteristiska ekvationen
det(ÅM - K)
=
Åm-3k k Åm-2k k 0
k
0
k Åm-3k
Egenvärdena är Å1 = k/m, Å2 = 3k/m och Å3 = 4k/m, och systemets egenfrekvenser blir
Motsvarande modvektorer blir
Här har vi valt modvektorerna så att de uppfyller ortonormalitetsrelationerna, men normeringen är inte nödvändig i denna typ av problem, eftersom amplitudfaktorerna ak står till fritt förfogande. Den allmänna lösningen
till det homogena systemet kan alltså mer explicit skrivas som
xi(t)
=
x2(t)
=
x3(t)
=
a1
a2
a3
r,;
~ cos(wot
+ 6i) + rn= cos( v 3wot + 62) + r..= cos(2wot + 63)
~ ~ cos(wot
+ 6i)
v6m
v6m a1
~ cos(wot
v6m
v2m
+ 61)
v3m
-
a2
r,;
~
r..= cos(2wot + 63) v3m a3
rn= cos( v 3wot + 62) + r..= cos(2wot + 63) v2m v3m
där wo = film. Vid principalsvängningama svänger alla tre massorna i takt och deformationerna är av formen
7.8. GEOMETBJSK TOLKNING AV KVADRATISKA EKVATIONER 143
7 .8
Geometrisk tolkning av kvadratiska ekvationer
En yta som kan representeras av en ekvation med formen ½xTAx-xTb+c=O
(A=AT/0)
kallas för en andragradsyta. Om x är tvådimensionell så talar man förstås om en andragradskurva. Vi skall här ge en översikt över de olika möjliga geometriska formerna för andragradsytor och andragradskurvor, genom att välja koordinatsystem så att ekvationen får en enkel form. Om A är inverterbar, så finns en entydig lösning till ekvationen
Ax= b Den entydigt bestämda punkten x kallas för ytans medelpunkt. Vi skall först ta fram de möjliga utseendena av en andragradsyta med medelpunkt. En flyttning av origo till medelpunkten, y = x - x, överför ekvationen i
½r T Ay + c' = o. Visa detta som övning och ange uttrycket för c'. Om c' = 0 så är ytan antingen en kon eller också består den av en enda punkt. Om c' # 0 så kan vi skriva ekvationen på formen yTKy
=1
KAPITEL 7. KVADRATISKA FORMER
144
om vi sätter K
= -(1/2c')A. Välj nu en ortogonal matris Q som diagonaliserar K,
I de nya koordinaterna 9, y
= Q9.
\ -2 A1Y1
så blir ekvationen
\ -2 \ -2 + A2Y2 + · · · + AnYn = 1·
Vi skall tolka detta geometriskt då n = 2 och n = 3. Låt först n = 2. Då är ekvationen
Om A1 < 0 och A2 < 0 så uppfyller inga punkter denna ekvation. Om A1 > 0 och A2 < 0 så får vi en hyperbel
- )2 - ( Y2- )2 = 1 Yl (1/v'Ti° 1/~ . Y2
Om A1 < 0 och A2 > 0 blir också kurvan en hyperbel. Om A1 > 0 och A2 > 0 så är kurvan en ellips med halvaxlar 1/v'Ti° och 1/,/X; och huvudaxlar längs koordinataxlarna.
- )2 + (1/,/X; - )2 -- 1. Yl Y2 (1/v'Ai
Yl □
Vi sammanfattar.
Sats 7.9 Antag att A är en symmetrisk inverterbar 2 x 2-matris. Ekvationen YT Ay
=l
har dd följande geometriska betydelse A negativt definit: ingen A positivt definit: ellips A indefinit : hyperbel Huwdaxlarna till kurvan är riktade längs egenvektorerna till A. I ellipsfallet är längderna av halva huvudaxlarna l / JX;, där Ak är egenvärden till A.
7.8. GEOMETRISK TOLKNING AV KVADRATISKA EKVATIONER 145
Vilken blir den geometriska tolkningen i specialfallet att egenvärdena är lika, .Å 1 = .Å2 > 0? För 3 x 3-matriser blir resultatet analogt. Den geometriska betydelsen avgörs även nu av egenvärdenas tecken, och det finns fyra olika fall , nämligen 0, 1, 2 eller 3 positiva.
Sats 7.10 Antag att A är en symmetrisk inverterbar 3 x 3-matris. Ekvationen YT Ay
=1
har dd följande geometriska betydelse A har O negativa egenvärden : ellipsoid A har 1 negativt egenvärde : enmantlad hyperboloid A har 2 negativa egenvärden: tvåmantlad hyperboloid A har 3 negativa egenvärden: ingen Huvudaxlarna till ytan är riktade längs egenvektorerna till A.
Om A är positivt definit, så betyder ekvationen YT Ay
=1
en ellipsoid. En sådan kan beskådas i figuren nedan.
Om A däremot är indefinit så erhålls olika typer av hyperboloider. Är två egenvärden positiva och ett negativt blir hyperboloiden enmantlad (sammanhängande), medan om ett egenvärde är positivt och två negativa blir den tvåmantlad (med två komponenter).
KAPITEL 7. KVADRATISKA FORMER
146
I fallet att något egenvärde är dubbelt. får vi ett helt plan av egenvektorer. Ytan blir då en rotationsyta med rotationsaxel vinkelrät mot detta plan.
Exempel 7.13 Ekvationen 3x~ + 2x~ + 3x~ - 2x1x2 - 2x2x3
=1
betyder, enligt exempel 7.4 en ellipsoid med halvaxlar 1, 1/ v'3 och 1/2, riktade längs vektorerna (1, 2, 1). (1, 0, -1) respektive (1, -1, 1) . D Lägg märke till att om man endast vill ta reda på av vilken typ en andragradsyta är, så räcker det att bestämma egenvärdenas tecken, och då kan man använda Gaus.selimination enligt Sylvesters tröghetslag. Vill man dessutom bestämma ytans storlek och form, måste egenvärdena beräknas, och för att få fram dess orientering i rummet behöver vi även ta fram egenvektorerna. Om matrisen A inte är inverterbar så finns ingen medelpunkt, och då kan i regel inte förstagradstermerna i ekvationen elimineras helt. Minst ett egenvärde till A är noll i detta fall, och efter ett diagonaliserande variabelbyte y = Qx får ekvationen ett utseende av formen ~(,X1Y:
+ A2Y~) - b~y1 - b;y2 - b;y3 + c = 0.
Om både ,\ 1 och ,\ 2 är skilda från noll. så kan koefficienterna 14 och b'2 elimineras genom kvadratkomplettering (flyttning av origo). Om den tredje koefficienten f'3 ej är noll, så kan konstanttermen också elimineras, och ekvationen får den nya formen
7.9. REFERENSER
147
där µ; = Åif~- Om µ 1 och µ 2 har samma tecken, så beskriver detta en elliptisk paraboloid, om de har motsatt tecken en hyperbolisk paraboloid (en sadelyta). Om~ = 0 så erhålls diverse typer av cylindrar, likaså om även ). 2 = 0. För figurer hänvisas till Persson & Böiers (1988), sid 23-24.
7 .9
Referenser
1. Spanne (1995a): mer om kvadratiska former och matristeori
2. Strang (1988): grundläggande lärobok i lineär algebra 3. Strang (1986): många tillämpningar av kvadratiska former 4. Uhlhom (1986): lärobok i mekanik
Inllsningsuppgifter: Kapitel 7 1 Definiera begreppet kvadratisk form och ange hur en sådan brukar beskrivas på matrisform. 2 Hur transformeras matrisen för en kvadratisk form vid ett lineärt koordinatbyte? Villren är skillnaden mellan denna transformationsformel och motsvarande formel vid lineära avbildningar? 3 Vad menas med diagonalisering av en kvadratisk form? 4 Vilken är den i princip enklaste metoden att diagonalisera en kvadratisk form? Hur genomför man denna metod med hjälp av matriser? 5 Vilken annan standardmetod för matri-
ser kan användas för räknemässig diagonalisering av kvadratiska former? 6 Kan Gausselimination användas för att enkelt bestämma egenvärden till en matris? 7 Beskriv ett samband mellan kvadratiska former och (lokala) optimeringsproblem för funktioner av flera variabler. 8 Vad menas med att en kvadratisk form är positivt definit? Vilka andra typer är definierade? 9 Hur avgör man enklast om en kvadratisk form är positivt definit? 10 Vilket samband finns mellan klassifikationen av en kvadratisk form och tecknen på egenvärdena för dess matris?
148
KAPITEL 7. KVADRATISKA FORMER
8 Numerisk lineär algebra
8.1
Preliminärt innehåll
Om detta kapitel någonsin kommer är tveksamt, men då blir innehållet kanske • Vektor- och matrisnormer • Matrisserier • Konditionstal • Potensmetoden för egenvärden. Invers iteration • Hessenbergmatriser och tridiagonala matriser. • QR-metoden • Kvadratiska former och minimeringsproblem • Rayleighs princip Nu får jag hänvisa till fortsättningskursema i matristeori (Spanne 1995a) och i numerisk analys. En utmärkt läro- och handbok i numerisk matrisräkning är Golub & Van Loan (1989).
149
150
KAPITEL 8. NUMERISK LINEÅR ALGEBRA
9 Insignal-utsignalmodeller
Vid sidan om tillståndsmodellerna används en helt annorlunda typ av matematiska modeller för system, nämligen insignal-utsignalmodellerna. Ett och samma tekniska system kan ofta beskrivas med båda typerna av modeller, men de har olika fördelar och nackdelar och belyser olika egenskaper hos systemet, och försvarar därför båda sin plats. Insignal-utsignalmodellerna är mer bundna till lineära system än tillståndsmodellerna, men kräver för övrigt mindre information om systemets detaljer för att kunna användas. I detta kapitel ger vi deras grundläggande egenskaper och definierar ett antal viktiga begrepp, linearitet, tidsinvarians, stabilitet och kausalitet.
9.1
~Black box~-modeller av system
De systemmodeller som vi hittills arbetat med, alltså tillståndsbeskrivningar, har varit av formen .-----.-.. D i - - - - - - - .
= {!x~=Cx+Dw. Ax+ Bw
Antalet tillståndsvariabler, alltså systemets ordning, kan vara mycket stort. Tillståndsmodeller är mycket användbara vid allmänna överläggningar och vid numeriska räkningar på dator, men om systemet är av högre ordning än tre, så är det nästan omöjligt att genomföra räkningarna för hand. Ofta är man inte intresserad av alla systemvariabler utan bara av någon enstaka, som beskriver ändamålet med systemet. Den som använder ett färdigt kommunikationssystem är inte betjänt av att veta allt om det otal strömmar, spänningar och 151
152
KAPITEL 9. INSIGNAL-UTSIGNALMODELLER
andra tillståndsvariabler som förekommer i hans apparater. För användaren räcker det för det mesta att känna till hur de signaler som kommer ut ur systemet förhåller sig till dem som sänts in. För sådana ändamål skall vi införa en enklare typ av systemmodeller, där man inte alls behöver räkna med systemens inre tillstånd. I dessa modeller representeras systemet enbart med en relation mellan en insignal och en utsignal. Systemet S anses vara bestämt , om man talat om vilka insignaler som är tillåtna 1 , och om man till varje tillåten insignal w kan bestämma motsvarande utsignal y. Grafiskt kan vi säga att man ersätter den komplicerade figuren ovan med en »låda» (amerikansk teknikerslang » black box»)
J.___
_w ..
S
.
y=Sw
och vad som sker inuti denna låda, beskrivs inte alls i systemmodellen. Lägg märke till att det är frågan om en matematisk modell av ett tekniskt system, och att man talar om alla tänkbara insignaler, inte bara den som verkligen uppträder för det konkreta systemet. Vidare låter man tiden sträcka sig över alla reella tal. från -oo till +oo . I praktiken har dynamiken för ett givet verkligt system en tidsskala som avgör för hur länge sedan -oo var. Eventuella jurister bland läsarna kan finna en mer formell definition av insignalutsignalrelationer på sid 182. Det finns system sådana att utsignalen y(t) vid en viss tidpunkt t endast beror av värdet på insignalen w( t) vid samma tidpunkt. Sådana system kallas för statiska. För de flesta system beror y(t) på större delar av förloppet av insignalen w, ibland hela, och systemet kallas för dynamislct. System och signaler kan vara mycket skiftande natur.
Exempel 9.1 (Telekommunikation)
s
y(t) utsignal: ljud ur hörlur
insignal: tal i mikrofon w(t) Här händer en massa på vägen, som man inte vill modellera exakt, och oftast inte heller kan. D kommer i fortsättningen att vara otillåtet slarviga och oftast utelämna precisa beskrivningar av vilka insignaler som är tillåtna för ett givet system. 1 Vi
9.1. »BLACK BOX»-MODELLER AV SYSTEM
153
Exempel 9.2 (Bil på ojämn väg) • insignal: vägens ojämnheter • utsignal: bilens skakningar
För att få en fullständig modell, måste vi ange bilens dynamik (antytt i figuren) och med vilken fart den framförs . D
Exempel 9.3 (Dator) • insignal: program och indata • utsignal: utskrifter D
Exempel 9.4 (Sveriges ekonomi) • insignal: diskontot • utsignal: valutareserv Detta exempel tycks vara lika aktuellt 2 varje gång jag skriver om den här boken. Man kan i detta fall fråga sig om systemet överhuvudtaget går att modellera matematiskt. I varje fall verkar det inte som om existerande modeller räcker till för att styra det. D
Det är naturligtvis möjligt att systemet har flera insignaler och utsignaler. För de typer av system och de problem som vi kommer att behandla i denna kurs innebär det inga nya svårigheter att utvidga till flera insignaler. Däremot innebär en sådan utvidgning stora problem då man försöker styra system. System är ofta utsatta, förutom för avsedda insignaler, för okontrollerbara störningar. Hur sådana skall modelleras kan vi inte gå in på här. Ofta beskrivs störningar med hjälp av stokastiska processer. 2 Detta skrevs första gången 1982. De formella styrvariablema har ändrats sedan dess, men slutsatsen om styrbarbeten verkar ännu inte ha blivit inaktuell.
KAPITEL 9. INSIGNAL-UTSIGNALMODELLER
154
De matematiska metoder för systemanalys som vi kommer att behandla här är inte enbart användbara på system där insignaler och utsignaler är funktioner av tiden. Det är också möjligt att arbeta med system . där in- och utsignalerna beror av en rumsvariabel eller ibland av mer abstrakta variabler. Huvudsaken är att dessa variabler kan beskrivas med reella tal eller med vektorer i R". I det första exemplet är in- och utsignalerna funktioner, definierade på ett intervall.
Exempel 9.5 (Balkböjning, statiskt) • insignal : belastning w(x) • utsignal : utböjning y(x)
utsignal: utböjning
y(x)
D
I nästa exempel är in- och utsignalerna fält, alltså funktioner definierade i det tredimensionella fysikaliska rummet.
Exempel 9.6 (Fältteori) En laddningstäthet p(x) i rymden ger upphov till en elektrisk potential U (x) (i stationärt tillstånd) • insignal: w(x)
= p(x)
• utsignal: y(x)
= U(x) D
I nästa exempel är insignalen definierad i det fysikaliska rummet och utsignalen är en funktion (svärtning) definierad i ett plan (bilden) .
Exempel 9.7 (Optiskt system) • insignal: verkligt föremål • utsignal: fotografisk bild D
9.2. HUR REALISTISKA ÅR INSIGNAL-UTSIGNALMODELLER?
155
I det avslutande exemplet kan både insignal och utsignal betraktas som funktioner av två reella variabler.
Exempel 9.8 (Bildbehandling) • insignal: inläst bild • utsignal: förbättrad bild Systemet S nedan tar bort rörelseoskärpan urinbilden:
s
H
0
Vi kommer i huvudsak att hämta exempel och beteckningar från fallet med tidsvariabla signaler. men samma matematiska begrepp och metoder återkommer i de övriga fallen. Det är bara den fysikaliska tolkningen som varierar.
9.2
Hur realistiska är insignal-utsignalmodeller?
För att få en uppfattning om begränsningarna i den nya typen av systemmodeller måste man tänka över när en insignal-utsignalrelation ger en tillfredställande bild av ett system. Vi vet ju från tillståndsmodellema att systemets utveckling beror inte bara av insignalen utan även av begynnelsetillståndet. När är inverkan av detta senare försumbart? För att besvara denna fråga kan vi använda resulteten i stabilitetsavsnittet. För ett stabilt system beror, eher ett insvängningsförlopp, under vilket en transient dör ut, rörelsen endast av insignalen. För instabila system dominerar däremot i regel de fria svängningarna, vilka bestäms av begynnelsevillkor. För att kunna behandlas som »svarta lådor» bör alltså system vara stabila. Vi kommer emellertid att ibland även att se på instabila system på detta sätt, ty de kan ingå som komponenter i stabila system, stabiliserade t ex genom återkoppling. Förutom stabilitet behövs även andra villkor. Det kan hända att inuti lådan finns komponenter, som inte är kopplade till ingången eller utgången. Detta kan leda till oönskade beteenden, som ej märks på den variabel, som vi valt att använda som utsignal eller som ej kan påverkas av insignalen. Villkor som används för att utesluta detta är styrbarhet, som grovt innebär att alla moder påverkas av insignalen, och obsen,erbarhet, som betyder att alla moder påverkar utsignalen. Dessa villkor undersöks närmare i reglertekniken. Se tex Åström (1976) .
KAPITEL 9. INSIGNAL-UTSIGNALMODELLER
156
I en systemmodell ingår strängt taget inte bara själva relationen y = Sw mellan insignal och utsignal, utan man måste också ange vilka insignaler som är tillåtna. En »last» med oändlig totalmassa på en balk måste nog uteslutas. Resonans och besläktade fenomen gör att somliga typer av insignaler inte kan behandlas i en del systemmodeller. Matematiskt innebär detta att operatorn S har en bestämd definitionsmängd, och att inte alla funktioner kan användas som insignaler. I speciella tillämpningar brukar man inte nämna detta explicit, men problemet dyker upp allt emellanåt.
9.3
Två typiska system
För att ge en mer konkret bild av insignal-utsignalrelationer inför vi två enkla modellsystem, ett elektriskt och ett mekaniskt.
Exempel 9.9 (Lågpassfilter) I det enkla RC-6.ltret i figuren ges sambandet
i
:,ellan insign{al ;:h uts:gnal e~igt exempel 1.2
~ •
dt y(t)
= RCx + RCw = x(t)
Vt
w
,
C
i+ .:
I_J_112_
Y 0
Lösning av differentialekvationen ger
y(t)
+
R
= e- 0
t
= 0 är oftast oväsentligt.
Att kärt barn har många namn gäller i högsta grad för stegfunktionen. ofta kallad Heavisides stegfunktion efter den engelske elektroingenjören Oliver Heaviside. Förutom den här använda har jag bl a funnit följande: t(t), b(t). u(t), H(t), U(t), l(t), l?(t), 77(t), x(t), Y(t), u(t), < t >o de flesta i kursböcker från LTH. En stegfunktion kan representera en likspänning som kopplas in vid tiden t = 0. För andra inkopplingstider får man använda förskjutna (translaterade) steg, 60 (t) = 9(t - a): 60 (t)
= 9(t -
a)
={
1
a
t
då t > a då t < a.
1 O
Andra viktiga signaler är pulser, speciellt rektangelpulser.
9.3 (Rektangelpuls) Enhets(rektangel)pulsen med bredd a definieras ge-
DEFINITION
nom annars. Värdet i diskontinuitetspunkterna är oftast oväsentligt.
Arean under enhetspulsen är
1
I
+ 00
_
00
På(t)dt
=aa = 1
oberoende av värdet på a > 0, vilket motiverar namnet enhetspuls. Med hjälp av stegfunktioner kan pulsen skrivas
t
9.5. STEGFUNKTIONER OCH PULSFUNKTIONER
En enhetspuls som börjar vid tiden t
( )= {! 0
PA t - r
=r
161
ges av
rr t < r.
Denna funktion har (som funktion av t) utseendet: k(t, r)
1 RC
t
r
Insignal-utsignalrelationen för lågpassfiltret kan nu skrivas y(t)
=
l
+oo
-oc
k(t, r)w(r) dr.
Beräkna som övning kil.(t, r) och visa att kil.(t, r)
➔
k(t, r) då~ ➔ 0.
□
KAPITEL 9. INSIGNAL-UTSIGNALMODELLER
168
lmpulssvar uppkommer på många olika områden och har många olika namn. De är lika viktiga för system med rumsvariabler. Exempel är viktfunktion (reglerteknik), Greens funktion. fundamental lösning eller kärnfunktion (partiella differentialekvationer. fältteori. kvantmekanik). point spread function (optik och bildbehandling) etc.
Exempel 9.13 (Balkböjning)
~(x,17~
(avläsningspunkt) x
y (belastningspunkt)
Vid balkböjningen i exempel 9.5 betyder k(x, y) utböjningen i punkten x förorsakad □ av en enhetspunktlast placerad i punkten y.
Exempel 9.14 (Fältteori) Enligt den elektriska fältteorin kan sambandet mellan laddningstäthet p och elektrisk potential U (från exempel 9.6) i vakuum skrivas som U(x)
= 4:eo
fff lx ~ Yl
p(YJ dyidY2dY3
(i SI-enheter). Här kan k(x, y)
1 _1 ~ = -41rfo 1X - y
uppfattas som potentialen i punkten cerad i 'fi. Relationen Sw(t)
=
(i vakuum)
x förorsakad
av en enhetspunktladdning pla□
1_:
00
k(t, T)w(T) dT
innebär en betydande principiell förenkling. Systemets egenskaper är nu sammanfattade i en enda funktion k(t, T) av två variabler, impulssvaret. Analogin mellan diskret och kontinuerligt blir i detta fall diskret Y=AX
kontinuerligt y=Sw
n
Yi
= Lai;x;
y(t)
=
l+oo k(t, T)w(T) dT -oo
j::;:l
lmpulssvaret för ett lineärt system är alltså helt analogt med matrisen för en lineär avbildning.
9.9. TIDSINVARIANS OCH TRANSLATIONSINVARIANS
9.9
169
Tidsinvarians och translationsinvarians
Allmänna lineära system är trots förenklingarna ovan ganska besvärliga att hantera analytiskt. Vi skall därför inskränka oss till en speciell typ av lineära modeller. Många system har den egenskapen att deras inre egenskaper, deras dynamik inte förändras med tiden. Härav följer att systemets svar på en insignal av en viss form alltid är detsamma (vid start från vila), dvs om insignalen fördröjes en viss tid T så blir utsignalen av samma form som innan, fast fördröjd tiden T. DEFINITION 9.6 Systemet S kallas tidsinvariant om följande villkor är uppfyllt för alla tillåtna insignaler w och alla tider T: Om insignalen w(t) ger upphov till utsignalen y(t) så ger den fördröjda insignalen w(t - r) upphov till den fördröjda utsignalen y(t - r). Speciellt krävs att om w(t) är en tillåten insignal så är w(t - a) en tillåten insignal för varje a.
Definitionen illustreras i figuren nedan. w(t)
y(t)
s
t
t
w(t - T)
y(t - T)
s t
T
T
Motsvarande egenskap för system med rumsvariabler kallas för translationsinvarians. Fysikaliskt innebär translationsinvarians att systemet har samma inneboende egenskaper i varje punkt i rummet, att det är homogent i rummet. Nästan alla systemmodeller som vi skall arbeta med är tidsinvarianta eller translationsinvarianta. Men det finns många praktiska system som inte uppfyller villkoren ovan.
Exempel 9.15 (Amplitud.modulation) Amplitudmodulation med undertryckt bärvåg är en lineär operation som kan beskrivas med formlerna y(t)
= a cos(w0 t)w(t)
där w0 är bärvågens vinkelfrekvens. Detta system är lineärt men ej tidsinvariant, ty w(t - r)
➔
a cos(w0 t)w(t - r)
KAPITEL 9. INSIGNAL-UTSIGNALMODELLER
170
medan
y(t - r)
= acos(w0 (t -
r))w(t - r).
Dessa uttryck stämmer i allmänhet inte överens om r inte är en heltalsmultipel av bärvågens period. w0 r /21r = heltal. Utsignalen är beroende av insignalens fas. D
Exempel 9.16 (Sampling) Sampling är en lineär men inte tidsinvariant operation. Om till exempel insignalen är periodisk med samplingintervallet som period, så kommer alla sampelvärden att bli lika. Vilket sampelvärde som avläses beror på insignalens fas relativt samplingsintervallen, och en fördröjning av insignalen leder till ett annat resultat än samma fördröjning av utsignalen. D Exempel 9.17 En balk som är upplagd på eller inspänd i stöd med fixt läge i rummet är inte translationsinvariant. En förflyttning av lasten men inte av stöden ändrar formen på utböjningen. Nedan syns utböjningsformema för en fritt upplagd balk. belastad med en enhetslast i två olika punkter. Det som syns i figuren är alltså impulssvar k(x, y) för två olika värden på y. k(x,yi)
k(x,112)
~ D
Om impulssvaret för ett lineärt system är känt, är det lätt att avgöra om systemet är tidsinvariant. lmpulssvaret k(t, r) vid tiden r måste i så fall ha samma form som impulssvaret k(t, 0) vid tiden Omen är fördröjt tiden r:
k(t, r)
= k(t -
r, O)
Sätt
h(t)
= k(t, 0)
Då är alltså
k(t, r)
= h(t -
r)
för alla t och r
och systemet bestäms helt av »funktionen» h(t). Omvänt gäller att om impulssvaret endast beror av skillnaden t- r så är systemet tidsinvariant. Visa detta som övning. Utsaga 9.2 Varje lineärt tidsinvariant system S har en insignal-utsignalrelation av formen
y(t) = Sw(t) =
l
+oc
-oc
h(t - r)w(r) dr
9.9. TIDSINVARIANS OCH TRANSLATIONSINVARIANS
171
där h( t) kallas för systemets impulssvar. Omvänt är varje system med en insignalutsignalrelation av denna typ lineärt och tidsinvariant.
Exempel 9.18 (Impulssvar för RC-fllter) Vårt enkla första ordningens RCfilter är givetvis tidsinvariant, vilket bland annat framgår av formeln k(t, r)
= ;ce-(t-T)/RC9(t -
r)
från exempel 9.12, där vi ser att impulssvaret bara beror av skillnaden t kan alltså karakteriseras med sitt impulssvar vid tiden r = 0, h(t)
T.
Det
= - 1-e-t!RC9(t) RC
som återges i figuren. 1 RC
h(t)
t
T
Fysikaliskt betyder h(t) utspänningen då filtret utsätts för en kort normerad spänD ningspuls vid tiden t = 0. I resten av kursen skall vi arbeta med metoder för att analysera lineära tidsinvarianta system. Detta är tillräckligt för att täcka en mycket stor klass av modeller. • De tillståndsmodeller, som vi arbetat med i den tidigare delen av kursen, alltså
i princip av typen
{!x~==
Ax+ Bw
Cx+Dw där A, B, C och Där konstanta matriser, leder, som vi skall se i avsnitt 9.12, till lineära, tidsinvarianta och kausala insignal-utsignalrelationer. • För många system finns en relation av formen y
00
ltl < 1 ltl > 1
~dw w
1. Om vi tar realdelen av formeln, så visar
cos(wt~ sin(w) dw
263
13.5. INVERSIONSFORMELN
har värdet 1r då ltl < 1 och O då It! > 1. En lärdom av detta är att integralen ovan
inte beror kontinuerligt på parametern t, fastän integranden cos(wt) sin(w)/w gör ~-
D
Ett annat exempel är den dubbla dämpade exponentialfunktionen från exempel 13.3.
Exempel 13.6 Från exempel 13.3 kan vi hämta transformparet
= e-it
f(t) {
-
2
= 1 + w2·
f(w)
Inversionsformeln ger i detta fall
r+oo ekalt
J_
00
2
1 + w2
d
(:=:_) = e-lt:. 21r
Verifiera som övning denna formel med residykalkyl.
-
D
Ur inversionsformeln följer att f = :F f bestämmer f (väsentligen) entydigt. Om F = :Ff så säger vi att f är inversa Fouriertransformen till F och skriver
f Inversionsformeln kan då skrivas :F- 1 F(t)
=
= :f"-lF_
1:
00
ekaltF(w)d (~).
Alla regler för Fouriertransformationen gäller med små modifikationer även för den inversa Fouriertransformationen. Genom inversionsformeln förklaras symmetrin mellan räknereglerna 3 och 4 samt den mellan 6 och 7. Enda skillnaderna är ett tecken samt faktorerna 21r. Den naturliga integrationsvariabeln på frekvenssidan är inte vinkelfrekvensen w utan frekvensen v = w /21r, i tekniska sammanhang mätt i Hertz. Använder man sig av frekvensen som variabel, vilket är vanligt såväl i elektroteknik som i ren matematik, så blir emellertid derivationsreglerna mindre enkla. Symmetrin mellan direkt och invers Fouriertransformation kan uttryckas i formelparet :Ff(w)
= F(w)
:FF(w)
= 21rf(-w)
Denna symmetri medför att en tabell över direkta Fouriertransformer efter teckenbyte och insättning av en faktor 21r också kan användas som en tabell för att bestämma inversa Fouriertransformer.
KAPITEL 13. FOURIERTRANSFORMATIONEN
264
Exempel 13. 7 Fouriertransfonnen av funktionen
t
J(t)
= (l + t2)2
skall bestämmas. I en tabell har man hittat transfonnparet
Funktionen till höger påminner onekligen om den givna, men beror av w. Sätt
.
G(w)=-41
w (1
+ w2 )
2-
lnversionsfonneln säger då att
Härav följer att D
Svårigheten vid denna typ av räkningar består bara i att hålla reda på beteckningarna. Metoden skall dock inte användas i onödan i sådana fall då en direkt användning av tabellen lättare leder till målet. Nästa exempel är »inversen» till exempel 13.1. I det skall vi bestämma transformen av en sincfunktion.
Exempel 13.8 Ur transfonnparet J(t)
= 9(t + 1) -_9(t -
1}
{
F(w) = J(w) = 2smw w
ser vi att transformen av
F(t)
sin t = 2-t-
är
21r/{-w)
= 21r(9(-w + 1) -
9(-w - 1)).
Men rektangelpulsen /(t) är uppenbarligen en jämn funktion {rita upp den) så att J(-w) = f(w), och vi får alltså transformparet sint -t-
F
~
1r(9(w + 1) - 9(w - 1)).
13.8. BEVIS AV FOURIERS INVERSIONSFORMEL
265
Fouriertransformen av en sincfunktion är alltså en rektangelfunktion. Lägg märke till att transformen av sincfunktionen inte är kontinuerlig, utan har språngdiskontinuiteter i w = ± 1. Härav följer att sincfunktionen inte kan tillhöra L 1 . Genom en direkt jämförelse med den harmoniska serien (med termer n - 1 ) kan man visa att integralen
är divergent. Däremot konvergerar ju
!
+oc
-00
sint
-dt=1r.
t
D
Vad betyder egentligen inversionsformeln? Vi återgår till härledningen ur Fourierserieutvecklingen för den upprepade funktionen fr- Denna innehåller frekvenskomponenter med vinkelfrekvenser kfl = kåw, k = 0, ±1, ±2, ... , där n = 21r/T. Då T --+ +oo så kommer ingående frekvenser att ligga allt tätare, med avståndet åw. I det operiodiska gränsfallet får vi en »kontinuerlig» summa, det vill säga en integral, av frekvenskomponenter. En operiodisk funktion innehåller i allmänhet komponenter för alla frekvenser, inte bara sådana som är heltalsmultipler av en grundfrekvens:
Detta resultat är grundläggande för svängnings- och vågfysiken, liksom för signalteorin. Formeln kan ibland förefalla något paradoxal. I exempel 13.2 så är funktionen /(t) = 0 för alla t < 0, men den kan i alla fall byggas upp ur harmoniska svängningar (av alla frekvenser), vilka har konstant amplitud för både positiva och negativa t. Dessa delsvängningar har dock sådana amplituder och faser att de genom »interferens» släcker ut varandra för alla negativa t.
13.6
Bevis av Fouriers inversionsformel
Vi inleder beviset med ett resultat av fristående intresse. Sats 13.3 Om
f tillhör L 1 , så gäller att
(a) J(w) är kontinuerlig och begränsad (av (b) (Riemauu-Lebesgue) J(w) --+ 0 då
11/111)
lwl --+ oo.
KAPITEL 13. FOURIERTRANSFORMATIONEN
266
Bevis: (a) Eftersom ie-i..lt /(t)I
= 1/(t)I
och e-i..it beror kontinuerligt på ;.,.J, så följer (a) direkt ur majorantsatsen. {b) Om 9(t) är en rektangelpuls så följer ur exempel 1 att g(w)-+ 0 då lwl -+ oo. Samma sak gäller för trappfunktioner, eftersom de är ändliga lineärkombinationer av rektangel pulser. Om f tillhör L 1 , så finns till varje€ > 0 en trappfunktion 9t sådan att
Vidare är för alla reella w
Välj a så att IBt(w)I < €/2 då lwl > a. (Detta är möjligt eftersom 9t är en trappfunktion, så att dess transform går mot noll i oändligheten.) Då är
om
lw I > a.
Alltså är
....
lim /(w)
lwl-+oo
= 0.
D
Vi skall nu ge ett bevis för inversionsformeln som gäller under förutsättningarna (a) / tillhör L1 (b) / är kontinuerlig i punkten t
= t 0 och har vänster- och högerderivata där.
För att förenkla räkningarna inför vi hjälpfunktionen
9(t)
= f(to + t) -
f(to)e-t 2 •
Denna tillhör L 1 samt är kontinuerlig och har vänster- och högerderivata i t = 0. Dessutom är 9(0) = 0. Vi inför ännu en hjälpfunktion 91 genom 91 (t) = 9(t)/t, så att 9(t) = tg1 (t). Funktionen 91 har både vänster- och högergränsvärde i t = 0, ty 91
(t) -_ 9(t) -_ 9(t) - 9(0)
t
t
13.8. BEVIS AV FOURIERS INVERSIONSFORMEL
267
och gränsvärdena är lika med motsvarande derivator för g it = 0. Alltså är g1 begränsad i en omgivning av t = 0, och utanför denna omgivning är den i L 1 , eftersom g(t) är det. Enligt frekvensderivationsregeln (7) följer då att g1 är deriverbar med -( gw ) = .ägt zdw och alltså ger insättningsformeln
Låt nu b--+ +oo och a--+ -oo. Ur Riemann-Lebesgue följer då att
r+oc }_
00
g(w) dw
= 0.
Vi skall nu återgå till den ursprungliga funktionen f. Den kan skrivas J(t
+ to) = g(t) + /(to)e-t2
och ur fördröjningsregeln (3), lineariteten och exempel 4 finner vi att
Alltså är
r~OC)
J_oc
eif-Ow J(w)
dw
r+oc
= J_oc
r+oc
g(w) dw + J(to)~ J_oo
e_..,2;4
dw
= 0 + /(to)~2~ = 21r/(to) vilket fullbordar beviset.
D
För funktioner, som inte är kontinuerliga i den betraktade punkten t = t 0 utan har en språngdiskontinuitet där, konvergerar inte inversionsintegralen absolut. Det är dock lätt att modifiera inversionsformeln så att den täcker även sådana fall. Vi ser först på ett enkelt exempel, en pulsmotsvarighet till en dipol.
KAPITEL 13. FOURIERTRANSFORMATIONEN
268
Exempel 13.9 (Dipolpuls) Låt / vara pulsen
/(!)={-i
f(t) 1--
då O < t < 1 då - 1 < t < 0 annars.
1
t
Dess transform beräknas lätt till
J(w)
= -4i sin2 (w/ 2)
om w =I= 0,
= 0.
J(O)
w
Lägg märke till att transformen J är en kontinuerlig begränsad udda funktion. som dock inte är integrerbar, eftersom den endast avtar som 1/ lw I för stora w. Funktionen / uppfyller villkoren i vårt bevis för inversionsformeln i alla punkter utom diskontinuitetspunkterna t 0 = 0, -1 och 1. För t 0 = 0 blir inversionsintegralen
och denna integral konvergerar ej som generaliserad integral i vår vanliga mening, eftersom delin_!egralerna från O till +oo och från -oo till O divergerar var för sig. Eftersom / är udda, är det dock rimligt att tänka sig att dessa båda oändliga delintegraler är lika stora och tar ut varandra. Man talar ofta om Cauchymedelvärden av integraler av denna typ, definierade (i enkla fall) genom
100
~00
CM _
F(w) (UiJ
=
0
~~ 00
10 F(w) _
0
(Ui/.
För en udda funktion är uppenbarligen Cauchymedelvärdet av integralen lika med noll, eftersom delintegralerna från O till a och från -a till Oexakt tar ut varandra. Vi konstaterar att för vår puls gäller en inversionsfonnel av typen
/(0 + 0); /(0 - 0) även i punkten t
= 0.
= CM 1:00 eiwO J(w)
d (~) D
Förfaringssättet i exemplet ovan fungerar generellt i punkter där / har en språngdiskontinuitet av enkelt slag, så att vänster- och högerderivatan existerar. Inversionsintegralen konvergerar då om den tas i Cauchymening. Dess värde blir inte säkert funktionsvärdet i språngpunkten, som ju kan väljas helt godtyckligt, utan medelvärdet av funktionens höger- och vänstergränsvärden.
13.7. FALTNINGSSATSEN
269
Sats 13.4 Antag att f tillhör L 1 • Om f har höger- och vänsterderivator i punkten t = t 0 sd gäller inversionsfonneln f(to
+ 0) + f(to -
0)
= CM 1..-oc eiwto l(w)
2
-oc
d ( ~). 21r
Bevis: Beviset av denna sats är nu mycket enkelt. Varje/ som uppfyller villkoren i satsen kan skrivas som summan av en dipolpuls som i exemplet ovan. translaterad till t = t 0 , och en funktion som är kontinuerlig i t = t 0 och för övrigt uppfyller villkoren i den tidigare inversionsformeln. För båda summanderna gäller formeln i satsen. o
f
Slutligen kan anmärkas att symmetrin mellan / och visar att inversionsformeln gäller punktvis även för vissa / som inte tillhör L 1 , till exempel /(t)
= -sint t
med
f(w)
= 1r (lJ(w + 1) -
lJ(w - 1)).
Lägg märke till att transformen av funktionen / är diskontinuerlig.
13. 7
Faltningssatsen
Redan vår introduktion av Fouriertransformationen visar att den hänger nära samman med faltningar. Själva definitionen kan ju skrivas som
Viktigt är att det finns ett mycket enkelt uttryck för Fouriertransformen av en faltning. Detta ges i
Sats 13.5 (Faltningssatsen) Om f och g tillhör L 1 sd är :F(/•g)
-
= :Ff. :Fg.
Bevis: Med användning av formeln ovan blir beviset en rent algebraisk räkning: f * g(w)
. . . . = e-twt ((f•g)•e""t) = e-awt (/•(g•etwt))
= e-twt
(/•(g(w)eiw')) = e-twt (/•eiwt) g(w) = f(w)g(w).
I det första steget har vi använt den associativa regeln för faltning. Denna gäller □ här, eftersom /(t) och g(t) tillhör L 1 , medan eiwt är begränsad. Insignal-utsignalrelationen för ett lineärt tidsinvariant system S
KAPITEL 13. FOURIERTRANSFORMATIONEN
270
__u·---·l._____s _
____,-y_=_S_w__
som i tidsområdet ges av den ibland komplicerade operationen faltning
övergår i frekvensområdet alltså i den mycket enklare operationen multiplikation y(w)
= h(w)w(w) = H(iw)w(w).
Beräkning av en faltning sker oftast enklast genom att man tar omvägen över frekvensområdet:
,
:,:-t
I
tidsområdet
:
-----~-----------------------------,----1
: :,:
frekvensområdet
----'---•I____:----· H(iw)·
_Y_.. H(iw)w
w
Fouriertransformering och invers Fouriertransformering sker, som vi redan nämnt, ofta lättast med hjälp av tabeller. Även numerisk beräkning av faltningar sker numera (sedan 1965) via en numerisk Fouriertransformation. Det finns nämligen en mycket effektiv algoritm för numerisk beräkning av Fouriertransformer, den s k snabba Fouriertransformen (FFT}. Denna finns kort beskriven i Spanne (1993) och något mera utförligt t ex i Strang (1986) eller Bracewell (1986). Exempel 13.10 Beräkna/•/, där /(t)
= e-t
2
_
Lösning: Vi använder först faltningsregeln:
/(t) ~ '✓'Ke-w2/4 F
-
/ •f(t) ~ /(w) 2
= 1re-w 12 2
13.8. PARSEVALS FORMEL
Eftersom e-1112 12 a = 1/,/2:
= e- 0. Ur skalningsregeln följer att
D
13.10. FOURIERTRANSFORMATION AV DELTAFUNKTIONER så att deltafunktionen är jämn. Detta bör vara intuitivt självklart. Sambandet mellan vinkelfrekvens w och vanlig frekvens II är ju w 11
= w /211'. Härav får vi formeln
279
= 211'11 eller
Transformen av en ren svängning kan alltså skrivas 1 ~ 6 (~)
eiwot ~ 6 ( w ;11'wo) Exempel 13.12 Beräkna Fouriertransformen av /(t)
= cosw0 t.
Lösning: Eulers formler ger
och alltså är □
Med hjälp av derivationsreglema, som fortfarande gäller i distributionsfallet, får vi ännu fler transformer. En dipol 6' ger :Fo'(w)
=
(derivationsregeln)
= iw:Fö(w) = iw · 1 = iw.
Upprepning av resonemanget ger regeln
Lineärkombinationer av delta i origo och dess derivator har alltså Fouriertransformer som är polynom i frekvensen. lnversionsformeln låter oss nu transformera även polynom it.
Exempel 13.13 Beräkna Fouriertransformen av funktionen /(t) Lösning: Vi kan skriva /(t)
= t • 1 och frekvensderivationsregeln
= t. ger □
280
KAPITEL 13. FOURIERTRANSFORMATIONEN
En upprepning av resonemanget i exemplet ger regeln
t" ~ i"2m5'"'(w) Sammanfattande ser vi att vid Fouriertransformation övergår lineärkombinationer av deltafunktionen i origo och dess derivator i polynom, och på grund av symmetrin i inversionsformeln gäller även omvändningen. Allmännare kan varje (lokalt integrerbar) funktion som växer högst som en potens av t för stora t Fouriertransformeras i distributionsmening. Om
1/(t)I :::; CltlN för stora t (N > 1) så är funktionen g( t) definierad av
g(t)
1
= 1 + t2N f(t)
absolut integrerbar (tillhör Li) och har alltså en kontinuerlig Fouriertransform i vanlig mening. Ur relationen
f(t) = (1 + t 2N) g(t) = g(t) + t 2Ng(t) och derivationsregeln på frekvenssidan följer att J(w)
= g(w) + (-l)N ~;Ng(w).
-
Denna formel bestämmer entydigt /(w) som distribution. Mängden av distributioner som har en Fouriertransform i distributionsmening (de så .kallade temperemde distributionerna) är inte uttömd med detta, men vi nöjer oss tills vidare med de ovan angivna exemplen. Se närmare Richards & Youn (1990) eller Friedlander (1982) för en utförligare framställning. Vi skall nu ännu en gång se på sambandet mellan Fourierserier och Fourierintegraler. Om f(t) är periodisk med period T så kan vi utveckla/ i Fourierserien
-00
Denna serie konvergerar åtminstone i distributionsmening, och det går att visa att den kan Fouriertransformeras termvis:
Detta kan i ord uttryckas så:
13.10. FOURIERTRANSFORMATION AV DELTAFUNKTIONER
281
Utsaga 13.1 Fouriertransfonnen av en periodisk funktion är ett tåg av likfonnigt utplacerade impulser, vars storlek är proportionell mot funktionens Fourierkoefficienter. Omvänt är inversa Fouriertransfonnen av ett sådant impulståg en periodisk distribution.
f(t)
f(w)
-:F
w
t
En funktion som är exakt periodisk med grundvinkelfrekvens n innehåller alltså skarpa frekvenskomponenter med vinkelfrekvenser kO, där k = 0, ±1. ±2, ... . I spektrum av naturligt uppträdande funktioner kan man ofta observera något liknande, men »verkliga» signaler är av begränsad varaktighet i en eller annan mening och det medför att »spektrallinjema» kring kO flyter ut och får positiv bredd. Vi skall nu se på Fouriertransformen av en ofta förekommande periodisk distribution.
Exempel 13.14 (Diracsta.ket) Beräkna Fouriertransformen av det periodiska impulståget
III(t)
i~
..-oo
IIl(t)
= L t5(t -
kT).
-oo
-T
T
2T
t
Denna distribution (ett » Diracstaket») är mycket användbar i tekniska sammanhang och kan väl vara värd en egen beteckning. Den här använda är stulen ur en amerikansk lärobok i bildbehandling. Varför den är lämplig förstår man kanske om man försöker rita upp 'grafen' av III(t). Egentligen är det fråga om en rysk bokstav ( med uttal lf al, på amerikanska shah) .
Lösning: Vi bestämmer först Fourierserien:
c1(III) = -1 T
1
e-ikOtIIl(t) dt = -1 per T
!T
e-ikOtIII(t) dt = -1 -T/ 2 T /2
1 e-ikOtt5(t) dt = -. -T/2 T
/ T.'2
Alltså är (eftersom varje periodisk distribution är summa av sin Fourierserie) 1
III(t)
=T
L eikOt •00
-oc
KAPITEL 13. FOURIERTRANSFORMATIONEN
282
och Fouriertransformation ger som ovan l ... oc ll(w) = T
+oo
L 2m5(w -
n
kO) =
-ex:
L 6(w - kO). -00
Fouriertransformen blir alltså ett impulståg i w-variabeln av precis samma typ, med D period n = 21r /T. Formlerna i exemplet ovan är så intressanta, att vi sammanfattar dem i en sats.
Sats 13.10 Antag att T > 0 och att TO +oo
L
= 21r.
l 6(t - kT) = T
k=-oc
k=-oc
+oo
L
eiknt
k=-oo -oo
+oo
L
Då är
6(t - kT) ~
n
L
6(w - kO)
k=-oo
Dessa formler har många användningar i samplingsteori och i andra sammanhang.
13.11
Två.dimensionell Fouriertransformation
Fouriertransformationen är nyttig inte bara för signaler som är funktioner av tiden, utan även för sådana som beror av rumsvariabler. Ett viktigt fall är bilder. En bild kan betraktas som en funktion av två variabler, där funktionsvärdet är svärtningsgraden i en punkt på bilden. En färgbild representeras i regel av ett vektorfält, där de tre vektorkomponenterna svarar mot t ex mättnaden av rött, grönt och blått eller av deras komplementfärger. Ett viktigt hjälpmedel vid bildbehandling är den tvådimensionella Fouriertransformationen, definierad av
Här ser vi bilder av författaren och hans Fouriertransform (egentligen bara absolutbeloppet):
13.12. HUR FOURIERTRANSFORMERAR MAN I PRAKTIKEN?
:F
283
>
För närmare information om användning av Fouriertransformationen i bildanalysen hänvisas till Gonzalez & ~•oods (1992) . Två viktiga användningsområden är bildrekonstruktion och bildkompression . Det kan också nämnas att de nu så kända magnetkamerorna. som använder kärnspinnresonans (~MR) för att fä bilder av det inre av människokroppen, bygger på tvådimensionell Fouriertransformation. Det som fysikaliskt mäts är Fouriertransformen av en funktion , väsentligen protontätheten i ett plant snitt av kroppen, och bilden fås genom en numerisk invers Fouriertransformation av mätdata. Denna bygger på inversionsformeln
och på den snabba Fouriertransformationen. Bland andra tillämpningar av två- och flerdimensionell Fouriertransformation nämner vi bara vågrörelseproblem inom mekanik. elektromagnetism och kvantmekanik och generellt lösning av partiella differentialekvationer. Den tvådimensionella Fouriertransformen har en direkt fysikalisk innebörd i optiken (Fraunhoferdiffraktion) som ibland har använts för analog beräkning av Fouriertransformer.
13.12
Hur Fouriertransformerar man i praktiken? De allmänna sambanden mellan en funktion f (t) och dess frekvensspektrum f (w) ges av Fouriers formler
-
KAPITEL 13. FOURIERTRANSFORMATIONEN
284
I praktiska situationer kan den ena funktionen vara given och man behöver bestämma den andra. Då finns ett antal olika metoder att välja på. Vilken man använder beror på hur funktionen är given och på vilka hjälpmedel man har till hands. Funktionen kan vara bestämd av en formel. den kan vara given av en tabell av mätvärden eller numeriskt beräknade värden eller också kan den vara representerad analogt, t ex som elektrisk spänning, optiskt eller magnetiskt. Metoderna kan indelas i två huvudklasser, en som bygger på matematisk beräkning av integralerna och en som bygger på fysikaliska effekter som kan beskrivas med Fouriertransformationen.
A. Beräkning av integralen 1. Insättningsformeln
Integranden har endast i extremt enkla fall en elementär primitiv funktion och vi har ovan angett de viktigaste fallen. 2. Komplex integration och residykalkyl
Komplexa metoder är mycket användbara och är det viktigaste hjälpmedlet för att bestämma nya transformer. 3. Tabellslagning En stor mängd Fouriertransformer finns tabellerade och i regel behöver man inte själv beräkna de transformer man skall använda, om de kan uttryckas med bekanta funktioner. Det krävs dock goda kunskaper om Fouriertransformer för att man effektivt skall kunna utnyttja tabellverken. 4. Numeriska metoder För numeriskt givna funktioner och då inga enkla uttryck finns för transformen så måste integralerna beräknas numeriskt. Det finns mycket effektiva numeriska algoritmer (den snabba Fouriertransformen, FFT) för dessa beräkningar. Numera finns till och med speciella integrerade kretsar för ändamålet.
B. Analoga metoder 1. Bandpassfiltrering Ett smalbandigt bandpassfilter med mittfrekvens w 0 släpper endast igenom Fourierkomponenter av / med frekve~r nära w 0 och utsignalen från filtret ger alltså en möjlighet att mäta /(w0 ). Genom att sända signalen genom ett batteri av filter med olika mittfrekvenser, eller vanligare, genom att successivt svepa filtrets mittfrekvens, kan man mäta spektrum av /. Många instrument för frekvensanalys bygger på denna princip. Ett icke kommersiellt sådant återfinner man i varje öra. Det är där fråga om mekaniska filter (i basilarmembranen), kombinerade med en hel del
13.13. STABil.A SYSTEM I TIDS- OCH FREKVENSOMRÅDET
285
avancerad signalbehandling i hjärnan. Även ett spektroskop bygger på samma princip. Ljus med olika frekvenser, dvs med olika färger. sänds ut i olika riktningar i rummet av ett prisma eller ett gitter. 2. Optisk Fouriertransformering
t.
I ..
t.
t, fraurmofcr paucm
hcld amplitudc
.::.. f
:Ff
fil
J±L.~, "
2
0
J
.
fldd
I
.,t,J
!
O
L••.
M1 d
- -2-
2
:F:F I
.mp1,,uc1r
2
Transform oti;cctplanc
o, foul planc
Iffla9C
planc:
En lins fungerar med god approximation som en (tvådimensionell) Fouriertransformator (genom Fraunhoferdiffraktion). Detta kan utnyttjas för att snabbt Fouriertransformera en funktion av två variabler, och likaså för optisk filtrering av bilder. För en mer ingående behandling hänvisar vi till moderna läroböcker i optik, t ex Hecht ( 1987), kapitel 11 och 14, i vilka man även kan finna många andra tillämpningar av Fouriertransformationen. Allmänt kan sägas att med datorernas utveckling har de numeriska Fouriermetodema fått ett allt större användningsområde. De ersätter numera ofta äldre analoga metoder. Likaså har Fouriermetoder blivit ett praktiskt beräkningsmedel på många områden där de tidigare endast var en del av den grundläggande teorin.
13.13
Stabila system i tids- och frekvensområdet
Vi sammanfattar till slut de olika beskrivningar vi fått fram för ett stabilt lineärt tidsinvariant system.
KAPITEL 13. FOURIERTRANSFORMATIONEN
286
__u·--·l...___s_
____.-y_S_w__ ..
Sambandet mellan insignal w och utsignal y är i tidsområdet y(t)
r+oo
= h•w(t) = J_
00
h(t - r)w(r) dr
och i frekvensområdet y(w)
= H(iw)w(w).
Impulssvaret h ges av h(t)
= YtS(t)
och frekvensfunktionen H (iw) av
lmpulssvar och frekvensfunktion kan rekonstrueras ur varandra genom integralformlerna H(iw) h(t)
r+oo
= :Fh(w) = }_
00
e-iwth(t)dt
= :F- 1 H(t) = f+oo eiwt H(iw)
J_oo
w)
d(2
1r
Mätning av den ena av funktionerna h eller H gör det möjligt att beräkna den andra numeriskt. Detta används i många moderna spektroskopiska metoder (en typ av Fouriertransformspektroskopi) där impulssvaret (h(t)) mäts direkt i stället för spektrum, som sedan beräknas numeriskt. I en annan typ av Fourierspektroskopi mäts en funktion av typen S11 (se beviset av Parsevals formel, sid 271, autokorrelation) och en invers Fouriertransform ger sedan IH(iw)l 2 , alltså i princip amplitudspektrum. Insignal-utsignalrelationen kan uttryckas i formeln y(t) ~ Energiinnehållet»
. H(iw)w(w) d ( w) = Jr+oo _ e1u1t 21r 00
i en signal w ges av
13.14. REFERENSER
287
och »energitätheten» i frekvensområdet är
I många tekniska sammanhang är det lättast att mäta beloppet lw(w)I, medan argumentet av w(w) är svårare att komma åt, och ibland mindre intressant. Sambandet mellan insignalens och utsignalens »energitätheter» ges direkt av
och utnyttjar alltså bara amplitudfunktionen A(w) = IH(iw)I. Om insignalens amplitudspektrum lw(w)I är konstant så blir alltså utsignalens amplitudspektrum jy(w)I direkt proportionellt mot A(w), oberoende av insignalens fas. Vanlig spektroskopi kan tolkas så att man sänder in bredbandigt brus i ett atomärt system (genom t ex uppvärmning). På plåten får man efter uppdelning i frekvenskomponenter genom ett gitter eller prisma intensiteter proportionella mot IA(w)l2.
13.14
Referenser
1. Bracewell (1986): behandlar ur tillämpningssynpunkt många delar av Fourie-
ranalysen och dess angränsande matematiska områden. 2. Dym & McKean (1972): trevlig och ej alltför svårläst matematisk bok om Fourieranalys.
3. Erdelyi, Magnus, Oberhettinger & Tricomi (1953-54): stort tabellverk över integraltransformer. 4. Hecht (1987): lärobok i optik. 5. Friedlander (1982): läsbar matematisk introduktion till distributionsteorin och Fourieranalysen. 6. Gonzalez & Woods (1992): lärobok i bildanalys.
7. Richards & Youn (1990): ännu en läsbar matematisk introduktion till distributionsteorin och Fourieranalysen. 8. Fourier (1822): pionjärverket i värmeledningens teori. där också Fourieranalysen började.
9. Hörmander (1990): första bandet av standardverket i lineära partiella differentialekvationer.
288
KAPITEL 13. FOURIERTRANSFORMATIONEN
Inläsningsuppgifter: Kapitel 13 1 Definiera Fouriertransformen av en funk-
tion. Ange några villkor på funktionen som garanterar existensen av transformen. 2 Skriv upp räknereglerna för Fouriertransformationen. 3 Härled derivationsreglerna för Fouriertransformationen. 4 Härled Fouriertransformen av funktionen e-t2
5 Skriv upp Fouriers inversionsformel. 6 Vad händer om en funktion Fouriertransformeras två gånger? 7 Härled Fouriers inversionsformel. Vilka förutsättningar behöver vi lägga på funktionen för att härledningen skall fungera? 8 Hur kan impulssvaret för ett stabilt lineärt tidsinvariant system erhållas, om frekvensfunktionen för systemet är känd? Vilken är anledningen till kravet på stabilitet?
9 Härled faltningssatsen för Fouriertransformationen. 10 Vilken är den största nyttan med faltningssatsen? 11 Hur kan enkelt faltningar beräknas, om man har tillgång till datorprogram för Fouriertransformation? 12 Givet är frekvensfunktionen för ett lineärt tidsinvariant system. Ange och härled sambandet mellan insignalens och utsignalens Fouriertransformer. 13 Frekvensfunktionen för ett system kan bestämmas, om man känner en enda, nästan godtycklig insignal och svaret på denna. Beskriv hur och ange villkoret på insignalen. 14 Skriv upp Parsevals formel för Fouriertransformationen och ge en fysikalisk tolkning av den. 15 Härled Parsevals formel. 16 Vilket allmänt utseende har Fouriertransformen av en periodisk funktion?
14 Laplacetransformationen
14.1
Definition av Laplacetransformationen
Fouriertransformationen ger bland många andra tillämpningar samband mellan systembeskrivningar i tidsområdet och i frekvensområdet. En nackdel är att Fouriermetoder endast kan användas på stabila system (eller gränsfall av sådana), för vilka funktionerna eMolt är tillåtna insignaler. För att kunna behandla även andra typer skall vi införa en nära besläktad transformation som ger den allmänna överföringsfunktionen H (s) på samma sätt som Fouriertransformationen ger frekvensfunktionen H (iw). Sambandet
mellan impuls.war h( t) och överföringsfunktion H ( s) leder oss till att göra följande definition:
Låt f(t), -oo < t < +oo, vara en funktion. Dess Laplacetransform C,f(s) definieras av formeln DEFINITION 14.1
C.J(s)
=
l
+oc
-oc
e-•tJ(t) dt.
Laplacetransformen av/ är definierad för de komplexa s, för vilka den generaliserade integralen konvergerar. Vilka är dessa? Om s = f3) parallell med imaginäro axeln sådan att Laplaceintegralen för f konvergeror absolut i hela det inre av strimlan men inte någonstans utanför strimlan. Anmärkning: Området för konvergens och för absolut konvergens av Laplaceintegraler kan i allra högsta grad skilja sig åt. Ett exempel är funktionen /(t)
= et sin(1re')8(t)
> 0 men konvergerar absolut endast i Re s > 1. Dessutom gäller att Laplacetransformen F(s) kan fortsättas analytiskt till en funktion som är analytisk i helas-planet. (Allt detta kan man bevisa genom att partialintegrera två gånger.) För de enkla funktioner som vi oftast arbetar med (väsentligen rationella och nära besläktade) sammanfaller dock det inre av strimlorna för konvergens, absolut konvergens och analyticitet. där integralen konvergerar i Re s
Exempel 14.1 Om f(t)
= 9(t) så är
14.1. DEFINITION AV LAPLACETRANSFORMATIONEN
om Res > 0. Om Res < 0 så divergerar Laplaceintegralen. 1 Alltså är C8(s) = -, Res > 0. Konvergensgränserna är alltså a s +oo.
291
= 0 och
;J = □
Lägg märke till att stegfunktionen 8( t) inte har någon Fouriertransform i vanlig mening. Den kan däremot Fouriertransformera.s i distributionsmening. ty den är en begränsad funktion och speciellt tempererad.
Exempel 14.2 Om f(t)
= 8(t) -
1 så är
då Res< 0.
D
Alltså gäller att
8( t)
C.
1 s
.........+ - ,
Re s > 0 (1
och
C.
8( t) - 1 .........+
1 s
- ,
Re s < 0
I figurerna är definitionsstrimlorna (här halvplan) skuggade och transformens pol (is= 0) utmärkt med ett kryss. De båda transformerna ges alltså av samma formel Cf(s) = s- 1 , men har olika definitionsstrimlor. Vi skall senare se mer precist hur läget av definitionsstrimlan återspeglar uppförandet hos f(t) då t----+ ±oe. Ytterligare ett exempel på samma fenomen finner vi i nästa exempel.
KAPITEL 14. LAPLACETRANSFORMATIONEN
292
Exempel 14.3 Betrakta funktionerna h(t)
(e-t - et)(9(t) - 1)
1
1 -1
-1
t
t
1
t
1
Beräkning ger
2 Cf1(s)= l-s 2
2
Ch(s)
Res> 1
= 1-
C,!J(s)
s2
2
= 1 -s 2
Res< -1
-1 0
-
s 1
D
(dämpning) e'0 '9( t)
Res> Reso. So
S -
På motsvarande sätt får vi
Exempel 14.7 Beräkna Laplacetransformen av /(t) Lösning: Vi använder räknereglerna: 9{t) - 1
C, i----,
1
= e' '(9(t) 0
1).
Res< 0
-
s 1
(dämpning) e'
0' (
D
Res< Res 0 .
9( t) - 1) S -
So
KAPITEL 14. LAPLACETRANSFORMATIONEN
294
Linearitetsregeln är också nyttig.
= sinw0 t 8(t).
Exempel 14.8 Beräkna Laplacetransformen av f(t) Lösning: Enligt Euler är sin wot
1
.
= 2i (e""'
01 -
.
e-""'0
t )
Exempel 14.6 och linearitet ger .
.C(sm w 0 t 8(t))
1 1 = --: . 2z s - zu.1
1 0
1
2i s - (-iwo)
då Res> 0.
14.3
D
Laplacetransformation och analytiska funktioner
Fouriertransformationen har en mer direkt fysikalisk tolkning än Laplacetransformationen, men den senare har flera räknemässiga fördelar. En liten sådan är att det är lättare att skriva s än iw. Den väsentligaste är dock att om Laplaceintegralen konvergerar i en strimla, så blir .Cf(s) en analytisk funktion i det inre av denna strimla. Vi kan då använda kraftfulla metoder från funktionsteorin för att behandla transformen. Det gäller följande viktiga sats:
Sats 14.3 (Analyticitet av Laplacetransformen) Antag att Laplacetransfonnen F(s) av funktionen f(t) existerar i strimlan a Reso-
som gör det möjligt att transformera alla funktioner av typ polynom gånger exponentialfunktion gånger stegfunktion. D Det finns naturligtvis en tvilling till den bevisade formeln, med utseendet e, t" e•o'(9( t ) - 1 ) i----+
(
k! S -
So
)" 1 ,
+
Res< Reso.
14.5. INVERSIONSFORMELN
14.5
297
Inversionsformeln
Eftersom varje Laplacetransform kan tolkas som en Fouriertransform. så leder Fouriers inversionsformel till en integralformel för den inversa Laplacetransformen.
Sats 14.7 (Inversionsformeln för Laplacetransformationen) Om Laplacetransformen .CJ(s) = F(s) av funktionen J(t) existerar i strimlan o < Res < 8 och om o < u < {3, så är J(t)
= ~ f"+ioc e'tF(s)ds. 21rt lu-ioc
Gränserna u ± ioc betecknar här att integrationsvägen är den räta linjen s = u + iw, där u är fixt och w varierar från -oc till +oc. Formeln gäller punktvis för varje t till exempel om f är kontinuerlig och om integralen i formeln konvergerar. men är med lämplig tolkning, tex i distributionsmening, sann i det allmänna fallet. Härledning: Vi vet att
Fouriers inversionsformel ger då
e-"tf(t)
. F(u+iw)dw. = -1 1-TOC eswt 2,r
-oc
Det komplexa variabelbytet s = u + i...;, ds = idw (u fixt) överför integrationsvägen -oo < w < +oo på linjen Res = u parallell med imaginära axeln. Överflyttning av faktorn e-ut leder till formeln
f(t)
= -21
7r
lu+ioc u-ioc
1
euteiA.Jt F(s)-:- ds l
vilket fullbordar härledningen. □ Formeln har den tolkningen att f (t) kan skri~ som en summa av svängningar e•t = euteiA.Jt med konstant dämpning -u. Eftersom integranden e•t F( s) är en analytisk funktion så kan man använda Cauchys integralsats för att deformera integrationsvägen från en rät linje till andra lämpliga kurvor. Det är ofta möjligt att använda residykalkyl för att beräkna inversionsintegralen. Nyttan därmed är tvåfaldig. Dels får vi effektiva metoder att bestämma inversa Laplacetransformer, dels leder metoden till fysikaliskt viktiga samband mellan tidsfunktionen J(t) och singulariteterna, till exempel polerna, av transformen F(s) i det komplexa frekvensplanet. Speciellt för partiella differentialekvationer gör detta funktionsteori till ett oumbärligt hjälpmedel. En av de viktigaste följderna av inversionsformeln är att Laplacetransformen F(s) (med angiven definitionsstrimla) bestämmer funktionen J(t) entydigt, åtminstone i kontinuitetspunkter. Detta har man nytta av även om man inte behärskar komplex integration.
298
KAPITEL 14. LAPLACETRANSFORMATIONEN
Sats 14.8 Antag att fi(t) och h(t) har Laplacetransformerna F1 (s) respektive F2(s), och att F1(s)
= F2(s)
på något intervall (av reella axeln), där båda är definierade. Då är
i alla (gemensamma) kontinuitetspunkter till
/i
och
h-
Om funktionen F(s) (med given definitionsstrimla) är en Laplacetransform, så kommer vi att använda beteckningen
f(t)
= c-• F(t)
för den (~ntligen) entydiga funktion f(t), som uppfyller F = Cf och tala om den som den inversa Laplacetransformen av F. Av den direkta Laplacetransformationen ärver den en mängd räkneregler, som fås genom att man läser formlerna för ~ från höger till vänster.
14.6
Funktionsteoretisk härledning av inversionsformeln
För att visa nyttan av funktionsteori skall vi nu ge ett alternativt bevis för inversionsformeln för Laplacetransformationen. Beviset bygger på Cauchys integralsats. Vi härleder regeln endast för relativt regulära funktioner. Som i tidigare bevis räcker det att bevisa inversionsformeln för t = 0, på grund av förskjutningsregeln. Dessutom kan vi anta att u =Res= 0, på grund av »moduleringsregeln». Vi förutsätter att funktionen / är två. gånger deriverbar samt att /, /' och /" är (absolut) integrerbara. Inför funktionerna F+(s)
= fo+oo e-•'J(t)dt,
F_(s)
={ 0
00
e-•'J(t)dt,
Res?: 0 Res~ 0
Båda dessa funktioner är kontinuerliga och begränsade i de angivna halvplanen, och dessutom är F + analytisk i Res > 0 och F _ i Res < 0. Partialintegration två. gånger ger formlerna
Alltså är
14.7. LAPLACETRANSFORMATION OCH DISTRIBUTIONER då
299
lal ➔ +oo i höger resp vänster halvplan. Vi beräknar nu den komplexa integralen:
Här har vi bytt integrationsvägar enligt Cauchys integralsats2 , från ett intervall på imaginära axeln till en halvcirkel C+ i höger halvplan respektive en halvcirkel C _ i vänster halvplan.
iR
-iR Låter vi nu R
➔
oo så får vi den önskade formeln +ioc / -ioc
14. 7
F(s) ds
= 21ri/(0).
Laplacetransformation och distributioner
Det är möjligt att definiera Laplacetransformationen för en hel del distributioner, fast inte alla. Vi skall här nöja oss med deltafunktionen och dess derivator. Alla de räkneregler som vi angett för Laplacetransformationen gäller fortfarande i distributionsfallet. Insättning i definitionen av transformationen ger genast reglerna 2 Skall man vara riktigt noga, så har vi använt en variant av Cauchys integralsats där man inte förutsätter att integranden är analytisk på själva integrationskurvan utan endast kontinuerlig.
KAPITEL 14. LAPLACETRANSFORMATIONEN
300
Sats 14.9 Laplacetransformen av en deltafunktion är en exponentialfunktion:
ö(t) ö(t - to)
H
1
H
e-to•.
Exempel 14.11 Beräkna Laplacetransformationen av
/(t)
= 3o"(t) -
26'(t - 1) + M(t - 2).
Lösning: Vi utgår från transformen av ö och tillämpar successivt räknereglerna: ö(t) ö1 ( t}
H
derivation fördröjning o'(t - 1) ö(t - 2) derivation ö" (t)
Svar: F(s)
= 3s2 -
2se-•
1
H
1=s e-•s e- 2• • 1 = e- 2•
H
S • S
H H
s·
= s2
+ 5e- 2'.
D
Vi kan också omedelbart bestämma den inversa Laplacetransformationen till ett godtyckligt polynom. Exempel 14.12 Bestäm den inversa Laplacetransformationen till F(s)
Svar: /(t)
= s 10 + 3s1 -
= .c- 1 F(t) = 0< 10>(t) + 3o( 7l(t) -
6s3
+ 11.
6o< 3 l(t) + llö(t).
D
Lägg märke till att om / (t) är en vanlig funktion och Laplaceintegralen konvergerar absolut, så är F(s) begränsad på varje linje Res= u inuti konvergensstrimlan. För distributioner kan F(s) = F(u + iw) däremot växa polynomiellt då lwl ➔ oo. Vi kan nu ge ett fullständigt svar på frågan om vilka funktioner som är Laplacetransformer (av en distribution). Sats 14.10 En funktion F(s) är Laplacetransformen av en distribution f(t) då och endast då det finns en strimla i s-planet, parallell med imaginära axeln, där F(s) är analytisk och uppfyller ett tillväxtvillkor F(s) = O(lslN) då lsl ➔ oo.
Att Laplacetransformen av en distribution som kan transformeras uppfyller villkoret kan vi inte visa här, eftersom vi inte genomfört den allmänna definitionen av (tempererad) distribution. (Det är det omtalade kontinuitetskravet som saknas.) Omvändningen är enklare och mer användbar. Vi bildar G(s)
1
= s N+ 2 F(s).
14.8. RATIONELLA TRANSFORMER
301
Inversionsintegralen för G( s) blir då absolut konvergent och definierar en funktion g(t), l lu+ioc g(t) = 21ft. u-ioc e•tG(s) ds, som uppfyller
Cg(s)
= G(s).
Sätt f(t)
= g I. Ur tabell finner vi svaret
.t:,~F(t)
= IS"(t) + IS(t) + ½et9(t) -
½e-t9(t).
D
14.8. RATIONELLA TRANSFORMER
305
14.8.2 Inverstransfonnation med residykalkyl Om F( s) = Q( s) / P( s) är rationell och ett äkta bråk, så kan inversionsintegralen
I (t) = -21 . 111-ioc e•t F (s) ds 11"l
11-ioo
lätt beräknas med residykalkyl . Fört < 0 så måste integrationskonturen kompletteras åt höger, för t > 0 åt vänster, ty le•tl = eR.e•t skall vara litet på den tillagda konturen. I figuren nedan är transformens konvergensstrimla markerad i grått.
lms
Ims
t0
X
X
Res
s
X
Observera den negativa orienteringen av integrationsvägen i fallet t < 0. Eftersom F(s) -+ 0 då !si -+ oc, så följer av Jordans lemma (i funktionsteorin) att integralen över den tillagda halvcirkeln försvinner i gränsen då radien får gå mot oändligheten, och residysatsen ger följande resultat .
Sats 14.11 Antag att f(t) har den rationella Laplacetransformen
Q(s)
F(s)
= P(s) .
Då är f(t)
=
-L
Res(e•tF(s)) , t < O
Re•>11
Formlerna är sanna även om F inte är ett äkta bråk, ty polynomdelen av F(s) ger ej bidrag till residyerna och ger i inverstransformen bara upphov till termer av typ o(t), o'(t), . .. , vilka ju försvinner då t =I= 0. Däremot bestämmer de i detta fall inte inverstransformen fullständigt . Residyberäkningen kan ske med metoderna från funktionsteorin: 1. Enkelpol s
= p:
(a) Res(e•tF(s)) •=p
= e1" ,_,.p lim(s -
p)F(s)
KAPITEL 14. LAPLACETRANSFORMATIONEN
306
(b) Res a=p
(e•' QP((s)) = e1" PQ(p() om P(p) = 0. S) I
2. Multipelpol s
p)
= p av ordning~ N.
(a)
Res •=p
(e•'
1 (s -
p )k
) tk-1 e1" - ( k - l) !
(b) Utveckla G(s) = (s - p)N F(s) kring s = p, G(s)
= Ao + Ai(s -
p)
+ · · · + AN-1(s -
p)N-I
+ 0(1s -
PIN)
Då är
varav följer formeln
Residyberäkningen är i detta fall helt ekvivalent med en partialbråksuppdelning av den rationella funktionen F(s) följd av en tabellslagning. l,; r formlerna ovan kan vi hämta en mängd information om inverstransformen /(t) utan att utföra några beräkningar. Antag att funktionen /(t) har den rationella Laplacetransformen F(s), och att denna är ett äkta bråk. Då gäller l.
f (t) är en summa av termer av typ polynom gånger exponentialfunktion gånger 9( t) eller 9( t) - 1, där varje pol till F( s) bidrar med en term.
2. En pol p av ordning v ger en exponentialfaktor ef1C och ett polynom av grad V -1. 3. Poler till höger om definitionsstrimlan ger bidrag för negativa t, alltså med en faktor 9( t) - l, sådana till vänster för positiva t, alltså med en faktor 9( t). 4. f(t) är kausal då och endast då inga poler till F(s) ligger till höger om definitionsstrimlan, alltså om denna är ett högerhalvplan. 5. /(t) är integrerbar då och endast då definitionsstrimlan innehåller den imaginära axeln i sitt inre.
14.8. RATIONELLA TRANSFORMER
307
6. /(t) är begränsad då och enda.st då den imaginära axeln ligger inuti eller på randen av definitionstrimlan, och alla poler på den imaginära axeln är enkla. 7. J(t) är kausal och integrerbar då och enda.st då definitionsstrimlan innehåller halvplanet Res 2: 0. Det händer ofta att man vet att en funktion är kausal, eller begränsad, och känner en formel för dess Laplacetransform. Resultaten ovan bestämmer då definitionsstrimlan.
Exempel 14.16 Bestäm den kausala, den begränsade och den antikausala ( = 0 då
t > 0) inversa Laplacetransformen till den rationella funktionen s
F( 8 )
= -(s_+_l_)_(s---1-)-(s---2-)
Lösning: F har enkla poler i s = - I , s = I och s = 2. Den kausala inverstransformen har en definitionsstrimla som ligger till höger om alla poler, den begränsade har en strimla som innehåller imaginära axeln och den antikausala ligger till vänster om alla per ler:
kausal begränsad antikausal
•.----
: kausal I
(J
Re s > 2 -1 < Res < I
Res < - I
Det finns i detta fall ytterligare en inverstransform, men denna har knappa.st någon naturlig användning. För beräkning av inverstransformema är det tillräckligt att bestämma residyema till e•'F(s) i alla poler. Eftersom dessa alla är enkla, så är det lätt gjort. Vi kan till exempel använda handpåläggning: Res•=-1(e•'F(s)) = lim._,,_ 1(s + I)e•tF(s) = -¼e-t . Res.=1(e•'F(s)) = -½et Res.= 2 (e•tF(s)) = }e2' De sökta inverstransformema fås nu ur regeln att poler till vänster om definitionsstrimlan ger bidrag fört > 0, medan poler till höger ger bidrag fört < 0. Detta
KAPITEL 14. LAPLACETRANSFORMATIONEN
308
ger
Jämför detta resultat med det som vi erhållit på mer elementär väg, utan användning av residykalkyl. i exempel 14.13. D Vi löser nu ett analogt exempel med en multipelpol. Här ger nog residymetoden något enklare räkningar. Exempel 14.17 Bestäm den begränsade funktion vars Laplacetransform är
F(s)
= (s -
s
I)(s + 2) 3
Lösning: F(s) har en enkel pol i s = l och en trippelpol i s = -2. Villkoret om begränsning ger att definitionsstrimlan skall innehålla imaginära axeln och den är alltså av formen -2 < Re s < l. Vi bestämmer residyema till e''F(s) med de angivna metoderna. Res,=1(e''F(s)) = lim,➔1(s - I)F(s) = f7 e'. ~ Bilda G(s) = (s + 2) 3 F(s) = s/(s - 1) och utveckla kring s = -2. Det enklaste sättet är att använda Taylors formel på G(s). Derivation och insättning ger 2 1 , G"( -2) = - 2 G(-2) = 3, G'( -2) = - 9 27
B
och alltså
G(s) Alltså är
2
= -3 2
F(s)
1 12 -(s + 2) - --(s + 2) 2 + O((s + 2) 3 ) 9 2 27 . 1
= 3(s + 2) 3
1
-
1
9(s + 2) 2
-
1 1 27 s + 2 + R(s)
där R(s) är analytisk is= -2. Residyreglema ovan ger resultatet ,, ( ~(e F s))
2t = ( 32 - 91 t 2
1 -2' 27 )e .
För den begränsade inverstransformen ger polen i s medan den i s = -2 endast ger bidrag för t > 0.
=l
bidrag endast fört < 0
14.9. LÖSNING AV DIFFERENTIALEKVATIONER
309
Sammanfattande får vi svaret
f(t)
= 271 et (8 ( t) -
1)
1 2 1 + (3 t - 9 t-
1 -2t 27 )e 8(t).
D
Alternativt kunde vi ha bestämt partialbråksuppdelningen för F(s) med någon av metoderna i envariabelkursen. Här används lämpligen ansatsmetoden, och efter en del räkningar fås uppdelningen 1
F(s)
1
= 27 s -
2
1
1
+ 3 (s + 2) 3
1
1
1
1
---9(s+2)2 27 s+2
Tabellslagning leder sedan till samma resultat som ovan.
14.9
Lösning av differentialekvationer
Laplacetransformationen leder till en enkel metod att lösa lineära differentialekvationer med konstanta koefficienter (och ibland med polynomkoefficienter). Enligt derivationsregeln övergår den relativt komplicerade operationen derivation i tidsområdet i den enkla operationen multiplikation meds i frekvensområdet. En differentialekvation
transformeras i ekvationen
a,,s"Y(s) + an_ 1 s"- 1 Y(s) + · · · + a1 sY(s) + aoY(s)
= F(s)
där Y(s) = .Cy(s) och F(s) = .Cf(s). Ur denna ekvation erhålls med rent algebraiska räkningar Laplacetransfonnen Y(s) av lösningen,
Y(s)
1
= ------F(s) ans"+ · · · + a1s + ao
och den sökta lösningen y( t) kan sedan bestämmas genom inverstransformation. Exempel 14.18 Finn den begränsade lösningen till differentialekvationen y< 3 l + y" - 2y
Lösning: Sätt Y
= b'(t),
-x,
< t < +oo.
= .Cy. Laplacetransformation ger s3 Y + s 2 Y - 2Y = .Cb' = s
och alltså
KAPITEL 14. LAPLACETRANSFORMATIONEN
310
..,,.;
)( Y(s)
=
1
5
s3
+ s2
-
2
1
-1
(T
)( -1
Den rationella funktionen Y ( s) har poler i s = 1 och s = -1 ± i . Eftersom y skall vara begränsad så blir definitionsstrimlan -1 < Res < 1. Inverstransformation med t ex residykalkyl ger efter förenkling svaret y(t)
= ~e'(8(t) ;)
1) + ~(- cos t + 3 sin t)e-'8(t) .
D
;)
Laplacemetoden är i själva verket inte helt obesläktad med den metod vi haft innan att lösa sådana ekvationer, nämligen passningsmetoden, men leder i regel till kortare räkningar. Det är viktigt att lägga märke till förskjutningsregeln även vid inverstransformation. En faktor e- 0 • på frekvenssidan svarar mot en förskjutning t := t - a på tidssidan. Exempel 14.19 Finn den kausala lösningen till ekvationen y
Lösning: Sätt Y
,,
+y=
{ 1, 0 < t < 2 0 annars.
= C.y. Då ger Laplacetransformation s2 Y
+Y
1 1 = C.(8(t) - 8(t - 2))(s) = - - e- 2•-.
s
s
Högerledet är definierat för alla s =/:- 0 (och singulariteten i s = 0 är häv bar) . Alltså är 1 1 Y(s) = (1 - e- 2•) - - - då Res> 0. s s2 + l Funktionen Y(s) är inte rationell, och partialbråksmetoden för inverstransformering måste modifieras (fast residymetoden leder direkt till rätt resultat) . Lösningen är dock enkel. Inför hjälpfunktionen F(s)
1
= -s
1 2
s +
1-
s
1
s
s2
+ 1·
Funktionen F har den kausala inverstransformen
f (t) = 8( t)
- cos t 8( t)
= (1 -
cos t )8( t) .
14.10. LINEÅRA SYSTEM OCH LAPLACETRANSFORMATIONEN
Men
.c- 1 (e-:u F)(t) = f(t y(t)
311
2) så vi får svaret
= (1- cost)9(t) -
(1- cos(t - 2))9(t - 2).
D
I exemplen ovan har vi haft differentialekvationer med en enda obekant funktion. Laplacemetoden lämpar sig dock särskilt väl för system av differentialekvationer, som den överför i lineära ekvationssystem beroende på en parameter s. Dessa kan sedan lösas med rent algebraiska metoder. I exempel 14.25 finns en genomräknad lösning av ett sådant system av differentialekvationer.
14.10
Lineära system och Laplacetransformationen
Faltningssatsen kan direkt överföras från Fouriertransformationen till Laplacetransformationen.
Sats 14.12 (Faltningssatsen för Laplacetransformationen) Om f, g och f•g kan Laplacetransformeras sd är .C(/•g)
= .Cf . .Cg.
Definitionsstrimlan för .C(/ •g) innehåller den gemensamma delen av definitionsstrimlorna för .Cf och .Cg.
Faltningssatsen har viktiga systemteoretiska konsekvenser. Vi har redan i inledningen fått ett första samband:
Sats 14.13 Överföringsfunktionen H(s) för ett lineärt tidsinvariant system S är lika med Laplacetransformen av impulssvaret h( t), H(s)
= .Ch(s).
Inversionsformeln visar nu hur impulssvaret kan rekonstrueras ur överföringsfunktionen. Insignal-utsignalrelationen för systemet kan ju skrivas som y
= h•w
och faltningssatsen leder till motsvarande relation i frekvensplanet.
Sats 14.14 Antag att ett lineärt tidsinvariant system har överföringsfunktionen H(s). Om insignalen w(t) och utsignalen y(t) har Laplacetransformerna W(s) respektive Y (s) sd är sambandet mellan dem Y(s)
= H(s)W(s).
Överföringsfunktionen kan alltsd f ds som kvoten mellan utsignalens och insignalens Laplacetransformer H( ) = Y(s) s W(s).
KAPITEL 14. LAPLACETRANSFORMATIONEN
312
Anm. I t ex reglertekniken använder man oftast »ensidig Laplacetransformation» (se avsnitt 14.11). Då måste man lägga till förutsättningen att systemet startar "i vila" vid tiden t = 0, för att formeln i satsen skall gälla. Start i vila medför att behövliga begynnelsevärden av y och w är lika med 0. Exempel 14.20 Systemet S är lineärt, tidsinvariant och kausalt. Om man sänder in insignalen w(t) = sin t 8(t) så kommer signalen y(t) = sin 2t 8(t) ut. Bestäm systemets im pulssvar och svaret på insignalen w ( t) = t sin t 8( t). Lösning: Eftersom W(s)
=
1 2
s +1
och
Y(s)
=
2 s +4 2
så blir överföringsfunktionen H ( s) = y (s) = 2 s2 + 1 = 2 - 3 2 . W ( s) s2 + 4 s2 + 4
lmpulssvaret är alltså
h(t) = C,~H(t) = 26(t) - 3 sin 2t 8(t). Om w(t)
= t sin t 8(t) så är W(s)
d (
= - ds
s2
1
)
2s
+ 1 = (s2 + 1)2'
och utsignalens Laplacetransfonn blir Y(s)
2s
s2
+1
= H(s)W(s) = (s 2 +1 )2 2 s2 +4
~~
= ( s2 + 1 s2 + 4) = ; ( s2 : 1 - s2 : 4) ·
Systemets svar på insignalen w(t) y( t)
= t sin t 8(t) blir alltså
= 34 (cos t -
cos 2t )8( t).
D
Anm. Den uppmärksamme läsaren bör ha reagerat mot några skenbara paradoxer i problemet ovan. Vi har ju tidigare visat att ur ett lineärt tidsinvariant system kommer endast ut frekvenser som redan finns i insignalen. Skenbart sänder vi (för t > 0) endast in vinkelfrekvensen w = 1 men får ut vinkelfrekvensen w = 2. Förldaringen är följande. Funktionen sin t 8(t) är inte en ren harmonisk svängning för alla
14.10. LINEÅRA SYSTEM OCH LAPLACETRANSFORMATIONEN
313
t, och dess Fourierutveckling innehåller i själva verket alla frekvenser. Systemets (stationära) svar på w(t) = sint, -oo < t < +oo bestäms av H(li) = 0 och är alltså y(t) = 0. Den utsignal y(t) = sin 2t 8(t) som vi observerar fört > 0 är en fri svängning, en egensvängning till systemet, som exciterats av inkopplingsförloppet i t = 0. Lägg också märke till att systemet inte är stabilt, eftersom H(s) har poler s = ±2i på imaginära axeln.
De allmänna egenskaperna hos Laplacetransformer ger en mängd information i systemfallet.
Sats 14.15 Om systemet S är kausalt så är dess överföringsfunktion H(s) definierad i ett högerhalvplan Re s > o ( eventuellt tomt). Om systemet S är stabilt så är dess över/öringsfunktion H (s) definierad, kontinuerlig och begränsad på imaginära axeln. Om systemet S är kausalt och stabilt. så är dess överföringsfunktion H ( s) analytisk då Res > 0 och kontinuerlig och begränsad då Res ~ 0.
Tyvärr kan man inte dra de omvända slutsatserna utan ytterligare villkor på systemet. I de två närmast följande exemplen beskrivs system av stor praktisk betydelse. där de omvända implikationerna inte gäller. Exempel 14.21 Låt ett system ha impulssvaret h(t) = e- 12 • Ett sådant system kallas för ett Gaussfilter (av speciell typ), och är mycket använt som modellsystem i signalteorin. Det är stabilt men inte kausalt. Överföringsfunktionen
är ändå analytisk i helas-planet (men ej begränsad i något högerhalvplan).
D
Exempel 14.22 I vårt värmeledningsproblem (i kapitel 12.6) var överföringsfunktionen av formen H(s) = e - ~ z (principalgrenen). Denna funktion är analytisk i Res > 0 och kontinuerlig även då Re s = 0, men den har en singularitet, en förgreningspunkt is= 0. D Det finns också system som endast kan analyseras i tidsområdet. Exempel 14.23 Vi ser på ett kausalt system med impulssvaret h(t) = e' 2 8(t). Utsignalen är väldefinierad för varje (t ex) kausal insignal genom faltningsformeln
314
KAPITEL 14. LAPLACETRANSFORMATIONEN
Systemet är naturligtvis instabilt, och dess impulssvar växer så snabbt att Laplaceintegralen inte konvergerar för något enda s. Inte ens distributionsteori kan rädda konvergensen. och slutsatsen blir att frekvensanalytiska metoder inte kan användas på detta system. D Om vi däremot inskränker oss till det elementära fallet med system av ändlig ordning, alltså med rationell överföringsfunktion, så blir allting mycket enklare.
Sats 14.16 Antag att S är ett lineärt tidinvariant system av ändlig ordning med den rationella överföringsfunktionen H ( s) = Q( s) / P( s). Då gäller att 1. S är kausalt då och endast då H ( s) är definierad i ett högerhalvplan Res > o.
2. S är stabilt då och endast då H(s) är definierad på den imaginära axeln Res = 0 och begränsad där. I reglertekniska sammanhang förutsätts alla system vara kausala och man behöver ett stabilitetskriterium. Detta ges då av
Sats 14.17 Antag att systemet S är kausalt och av ändlig ordning. Då är S stabilt då och endast då dess överföringsfunktion H ( s) = Q( s) / P( s) uppfyller följande villkor: a) gradP 2:: gradQ
b) alla poler till H ligger i vänster halvplan Res < 0. I nästa kapitel kommer vi att använda detta kriterium för att få fram praktiska metoder att avgöra stabilitet av system.
14.11
Den ensidiga Laplacetransformationen
den ensidiga Laplacetransformationen ... Den lytta, halta, skumögda krympling som heter L-transformmetoden, är inte bara oduglig vid de tvungna svängningarna, dvs. i verklighetens värld, utan den saknar överhuvud taget existensberättigande ... E. Hallen {legendarisk professor i teoretisk elektroteknik vid KTH}.
Om f(t) är en kausal funktion så kan Laplacetransformen av/ skrivas som
Kausala och kausalt avskurna funktioner är så vanliga att detta uttryck har en egen benämning.
14.11. DEN ENSIDIGA LAPLACETRANSFORMATIONEN
315
Om J(t), -oo < t < +oc är en funktion, så definieras den ensidiga Laplacetransformen C, 1 f (s) av f genom formeln DEFINITION 14.2
Cif(s)
~ C,(6/)(s) = fo
00
e_.,J(t)dt.
Definitionsmängden för den ensidiga Laplacetransformen är alltså alltid ett högerhalvplan Res > o (eventuellt tomt). De flesta Laplacetransformer. som man träffar på i litteraturen, är i själva verket ensidiga transformer, och vi måste därför ta upp deras speciella egenskaper. Det är framför allt två räkneregler som skiljer den ensidiga transformen från den dubbelsidiga, nämligen derivations- och faltningsreglema.
Sats 14.18 (Derivationsregeln för C, 1 ) C.1(/(nl)(s)
= s"C,if(s) -
s"- 1/(0)- · · · - sJ'"- 2l(0) - /(n-ll(0).
Härledning: Ur derivationsregeln för en avskuren funktion 6/(n)
= (6/)(n) _
/{0)tS(n-1) _ ... _ /(n-2i(0)t5' _ Jtn-ll(0)t5
får vi omedelbart regeln ovan genom tvåsidig Laplacetransformation, ty
C,(6/("l)
= C1(/1"l),
C,((6/)(nJ)
= s"C,(6!) = s"C,if
och
C,(t5(i:J)
= sk. □
Den långa svansen i formeln efter huvudtermen s" C, 1f uppkommer alltså därför att det språng som vi får då f skärs av i t = 0 måste deriveras n gånger. Formeln kan för övrigt också härledas utan användning av deltafunktioner, genom partialintegration n gånger. Den andra regel som är speciell för den ensidiga transformen följer omedelbart ur regeln för faltning av kausalt avskurna funktioner, (6!)•(6g)(t)
= 6(t)
(l' J(t - r)g(r) dr) .
Detta ger
Sats 14.19 (Den ensidiga faltningssatsen) C1
(l' J(t - r)g(r) dr) (s) = C,if(s)C1g(s).
Den ensidiga Laplacetransformationen kan även definieras för vissa distributioner. En förutsättning är att de är så pass regulära nära tiden t = 0 att produkten med stegfunktionen är definierad. Om f innehåller t5 och derivator av t5 belägna i t = 0. så uppstår det problem. Man måste då konventionsmässigt bestämma om dessa singulariteter skall ingå i produkten 6f eller ej. I elektroteknisk litteratur säger man då ofta att den undre integrationsgränsen är 0- resp 0+. Dessa bekymmer undvikes lättast om man konsekvent arbetar med den tvåsidiga transformen och distributionsderivator.
316
KAPITEL 14. LAPLACETRANSFORMATIONEN
14.12
Begynnelsevärden och ensidiga problem
. . . så är hans metod ( angiven tidigare av andra) att först L-transformera hela differentialekvationen, ett våldsamt ingrepp i denna. Man klipper av och kastar med berott mod bort hela det föregående förloppet eller tillståndet och får ersätta denna kunskapsförlust med konstanter, ja med förhoppningen att små tomtar komma flygande med begynnelsekonstanter i tillräckligt antal. E. Hallen.
En mycket vanlig typ av problem är sådana där man är obekant med eller ointresserad av förloppen innan en viss tidpunkt, som konventionellt kan väljas till t = 0. Informationen om det tidigare skeendet kan sammanfattas i tillstdndet vid tiden t = 0. vilket för differentialekvationer av ordning n kan uttryckas med n stycken begynnelsevärden. Man skaffar sig då en ekvation, som gäller för alla t, genom att multiplicera med en stegfunktion 9(t). Sedan kan denna ekvation behandlas med olika metoder. Den ensidiga Laplacetransformationen erbjuder färdiga metoder, ty den innebär just först avskärning av förloppet för t < 0 och sedan Laplacetransformation.
Exempel 14.24 Lös begynnelsevärdesproblemet
{
y" + 2y' + y = 4tet, t > 0 y(0) = -3, y'(0) = 2.
Lösning: Differentialekvationen gäller endast för t > 0, men efter multiplikation med 9(t) får vi en ekvation, som gäller även fört < 0, ty då är bägge leden lika med noll: 9y" + 29y' + 9y = 4tet9(t), -oo < t < +oo. Vi Laplacetransformerar denna ekvation. Sätt Y(s) regeln för den ensidiga transformen ger
C(9y") C(9y')
= C1y(s) = C(9y). Derivations-
= s 2 Y - sy(0) - y'(0) = s2 Y + 3s = sY - y(0) = sY + 3.
2
Alltså är
s2 Y + 3s - 2 + 2(sY + 3) + Y = C1(4tet)(s) = ( Ur denna ekvation kan vi lätt lösa ut Y:
(s 2 + 2s + l)Y(s)
= -3s - 4 + \
1 ) s-1 2
4
s-1
)2 .
14.12. BEGYNNELSEVÅRDEN OCH ENSIDIGA PROBLEM och alltså
Y(s)
-3s - 4 (s + 1)2
317
4 (s - 1) (s + 1)
= - - - + - - -2 - - -2
Vi kan nu finna fJy (men inte y) genom kausal inverstransformation. Partialbråksuppdelning ger 1 1 2 y ( 8 ) = - s - 1 + (s - 1) 2 - s + 1 och alltså är fJ(t)y(t)
= -fJ(t)e' + tfJ(t)e' -
28(t)e-'.
Men vi kan bara bestämma y fört > 0. Eftersom fJ(t) = 1 då t > 0 så får vi Svar: y(t) = -e' + te' - 2e-t, t > 0. Kontrollera alltid att begynnelsevärdena stämmer! Insättning i uttrycket för y ovan och dess derivata ger y(O)
= -1 + 0 -
2 = -3,
y'(O)
= -1 + 1 + 0 + 2 = 2. □
Stämmer!
Laplacemetoden är speciellt användbar på system av differentialekvationer. Vi illustrerar detta på ett exempel ur mekaniken.
Exempel 14.25 Två massor är förbundna med varandra och med fixa stöd genom lineära fjädrar, enligt schemat i figuren.
3m Yl
Y2
Fria svängningar av systemet beskrivs då av följande system av differentialekvationer: 3my~ = -15ky1 - 6k(y1 - Y2) { 2my; = -6kY2 + 6k(y1 - Y2)Vid tidpunkten t = 0 är systemet i vila. Vi ger då den vänstra massan en stöt så att den erhåller hastigheten y~ = 1 åt höger. Bestäm massornas vidare rörelse. Lösning: Efter hyfsning får vi systemet
KAPITEL 14. LAPLACETRANSFORMATIONEN
318
då t > 0. med begynnelsevärden y1 (0) = 112(0) = 0 (start i jämviktsläget) och y;(0) = 1. y; = 0. Sätt Y1 = C1Y1 = C(8y1) och Y2 = C.1112 = C(8y2)- Efter Laplacetransformation och insättning av begynnelsevärdena övergår differentialekvationerna i
Hyfsning av systemet ger k k (82 + 7 - ) Y1 - 2- Y2 = 1 m m
{
k
k
+ (8 2 + 6-)Y2 =
-3-Yi m
m
0.
Detta är ett lineärt ekvationssystem för Yi och Y2 , med koefficienter som beror av parametern 8. Det förmodligen enklaste sättet att lösa lineära ekvationssystem av denna typ och storlek är med användning av Cramers regel. Determinanten för systemet är 82
A(8)
=
k m
+7k
-3m
k m
-282
k m
k m
k m
k m
k m
= 84 + 13-82 + 36(-) 2 = (8 2 + 4-)(8 2 + 9-).
+6-
Enligt Cramer blir lösningen till systemet
Y1
I
= -A
I
Y2=A
1 0
k m
-282
k m
+6k m
82 +7- I k -3m
Det är lämpligt att införa beteckningen w
{
I 2 k = -(8 +6-) A m
0
= _!_3~ A m·
= ~- Partialbråksuppdelning ger
Y1=~ I 5 82 + 4w2
+~
5 82
I
+ 9w2
Y:-~ I -~ I 2 - 5 82 + 4w2 5 82 + 9w2
14.12. BEGYNNELSEVÅRDEN OCH ENSIDIGA PROBLEM
319
och invers Laplacetransformation i tabell ger sedan svaret
Y1(t)
= ~(sin(2wt) + sin(3wt)) aw
112(t)
= 3~
{
(9sin(2wt) - 4sin(3wt))
då t > 0.
D
Då man löst ett problem med Laplacetransformation så måste man kontrollera att resultatet uppfyller de givna begynnelsevillkoren. Om högerleden inte innehåller några deltafunktioner så kan detta göras direkt genom insättning i t = 0. I detta sammanhang måste påpekas att vissa begynnelsevärdesproblem för system av differentialekvationer saknar lösning, men Laplacemetoden ger till synes ändå ett svar, fast det är felaktigt. Detta svar är en lösning till differentialekvationerna, men har fel begynnelsevärden.
Exempel 14.26 Finn alla lösningar till systemet
+ y' + y = 0 x' + y' + 2y = 0 x'
{
x(0) = y(0) = 1.
"Lösning:" Sätt X(s)
= C.(9x), {
Y(s)
= C,9y.
Laplacetransformering ger
sX - 1 + sY - 1 + Y sX - 1 + sY - 1 + 2Y
=0 = 0,
alltså
{ Här är
sX + (s + l)Y sX + (s + 2)Y
=2
=2
a=lss s+ll=s s+ 2
och enligt Cramer är
X(s)
= ! 12 s
{ Y(s) =
s+1
2 s+2
!s I ss 21 = 2
Inverstransformering leder till
{
9(t)x(t) 9(t)y(t)
= 29(t) =0
I = ~s 0.
320
KAPITEL 14. LAPLACETRANSFORMATIONEN
och alltså {
x(t)
=2
y(t)
=0
för t ~ 0. Men detta är inte en lösning till det uppställda problemet, ty x(O) = 2 och y(0) = 0 uppfyller ej begynnelsevillkoren. I själva verket finns ingen lösning, ty ur de båda ursprungliga ekvationerna följer att y = 0, och det enda möjliga begynnelsevillkoret är y(0) = 0. Var ligger felet i lösningen? Jo, då vi sätter X(s) = C(8x) så antar vi att det D finns en lösning (x(t), y(t)), vilket inte är fallet i detta exempel. Exemplet är avsett att visa hur farligt det kan vara att blint använda en metod, om man inte har tillgång till precisa matematiska resultat, som bevisar att metoden säkert leder till ett korrekt resultat. Det visar också vad man kan vinna på en rimlighetskontroll av erhållna resultat. För lineära tillståndsekvationer med konstanta koefficienter kan det obehagliga fenomenet ovan inte inträffa (eftersom begynnelsevärdesproblem för sådana alltid har entydig lösning, vilket vi bevisat i kapitel 5) och de kan alltså riskfritt lösas med Laplacemetoder. Härvid Laplacetransformerar man ofta systemet skrivet på matrisform. Laplacetransformen av en tidsberoende matris erhåller man helt enkelt genom att transformera varje element för sig. Det mest användbara fallet är den ensidiga transformen av en exponentialmatris.
Sats 14.20 Den ensidiga Laplacetransformen av exponentialmatrisen eAi är resolventmatrisen RA(s) = (si - A)- 1 :
Konvergensområdet är halvplanet Res> u(A).
Bevisa denna regel som övning. Utgå tex från att
och använd derivationsregeln. Detta samband mellan resolventmatriser och exponentialmatriser kan vid handräkning ofta användas för att bekvämt beräkna exponentialmatriser.
Exempel 14.27 Beräkna e'A, där A är matrisen 0 -1 A= [ 1 0
J.
14.13. FALTNINGSEKVATIONER
321
Lösning: Vi beräknar först resolventen och inverstransformerar denna sedan för att få fram exponentialmatrisen.
RA (s) = ( s I -
A)- 1
s = [ - l
l i-i 1 [ 1l = [ 1+ 1 1] . = ~
8
ls -
8
52
8
s2
Kausal inverstransformering med tabell ger etA 9(t)
= .c,-iR
(t)
= [ c~t 9(t) sm t 9( t)
A
Alltså är etA
= [ c~s t
smt
- sin t
cost
-sint 9(t) cos t 9( t)
5 2-: 5
+1
s2 + 1
l
l
för positiva t. Men etA är en hel analytisk funktion av t, så att det följer ur identitetssatsen att denna likhet gäller för alla (även komplexa) t. eftersom den gäller för D reella positiva t.
14.13
Faltningsekvationer
Elementära modeller för tekniska system brukar innehålla differentialekvationer, men i många fall får man i stället integralekvationer. Vi skall inte här gå in på allmännare sådana ekvationer, men integralekvationer av faltningstyp kan ofta lösas med transformmetoder.
Exempel 14.28 Lös integralekvationen y(t) -
fot cos(t -
r)y(r) dr= cos t,
t
~ 0.
Lösning: Multiplikation med 9(t) ger en ekvation som gäller för alla t: 9y(t) - 9(t)
fot cos(t -
r)y(r) dr= cos t 9(t).
Den andra termen i vänsterledet känner vi igen som en (kausal) faltning av cos t 9(t) och 9y( t), och ekvationen kan skrivas
9y - (9 cos t) •9y Sätt Y
= 9 cos t.
= C( 9y). Laplacetransformation ger då Y - .C,((9cost)•9y)
= .C,(9cost).
KAPITEL 14. LAPLACETRANSFORMATIONEN
322
Enligt tabell och faltningsformel får vi
Y-
s
s2
+1
Y=
s
s2
+1
Här kan vi lösa ut Y och får Y(s)
(s - l) + !
s
=- - - = (s - 2)2 i2 23· s2 - s + I +i
Kausal inverstransformation ger nu att
Svaret blir alltså y(t)
v'a) '
v'3 + -1s m . -t = e' 12 ( cos-t 2 v'J 2
t > 0.
D
I sannolikhetsteorin och i reglertekniken (speciellt återkoppling) är integralekvationer med ensidiga faltningar vanliga. De kan ofta l ~ med hjälp av den ensidiga Laplacetransformationen (trots att de beskriver verkligheten, tex telefonväxlar).
14.14
Begynnelse- och slutvärdessatserna
En fullständig invers Laplacetransformation av även så enkla funktioner som de rationella kräver en hel del räknearbete, speciellt om det finns multipla poler och om nämnarens gradtal är stort. Det är därför av stort praktiskt intresse att kunna hämta information om den inversa transformen utan att behöva bestämma den fullständigt. I det kausala fallet kan man hämta sådan information ur begynnelseoch slutvärdessatsema. Vi formulerar dessa satser endast för rationella transformer, men de gäller med lämpliga, fast ofta invecklade, förutsättningar även i mycket allmännare fall.
= C,f är rationell. Om F är ett äkta bnlk, alltstl om nämnarens gradtal är större än täljarens, stl är lim J(t) = lim sF(s).
Sats 14.21 (Begynnelsevärdessatsen) Antag att f är kausal och att F
t ➔ +O
1 ➔ +00
Sats 14.22 (Slutvärdessatsen) Antag att f är kausal och att F nell. Om alla poler till sF( s) har negativ realdel, stl är
lim J(t)
t ➔ +oo
=
lim sF(s).
1 ➔ +0
= C,f
är ratio-
14.14. BEGYNNELSE- OCH SLUTVÄRDESSATSERNA
323
Bevis: Satserna gäller för transformparen
e1"9(t) Ä
tk-1
1 (s - p)"
(k - 1)!
och
Ä
~ 0. w
speglad Nyquist-._,,
Q
.......__
... \
I
w
•
= H(s)
•
w=O
Nyquist
Den beskrivs ofta i polära koordinater,
H(iw)
= A(w)ei,p(w)
med amplitudfunktion A(w) och fasfunktion 0.
15.S. NYQVISTDIAGRAM
335
Sats 15.5 ( N yquistkriteriet) Antalet instabila poler till det återkopplade systemet är lika med antalet instabila poler till det ursprungliga systemet plus det antal van, som Nyquistkun,an och dess spegelbild går runt punkten w = - I/ k i positiv led.
Vi skall här se på ett exempel på en analytiskt given överföringsfunktion och rita upp motsvarande Nyqvistdiagram. När diagrammet endast skall användas till varvräkning kring punkter på reella axeln, så är det tillräckligt att bestämma alla axelskärningar för H(C). Dessa fås ur nollställen till ReH(iw) och ImH(iw) samt gränsvärdet då lwl ➔ +oo.
Exempel 15.1 (Fasskiftsoscillatom) Filtret i figuren har överföringsfunktionen
C
C
C
rl
H s (RCs)3 ( ) - (RCs) 3 + 6(RCs) 2 + 5RCs + I
w
+ R
R
R
y
(övning i kretsteori). Det återkopplas negativt via en förstärkare med förstärkning k. För vilka k blir den återkopplade kretsen instabil? För sådana k kan den (kompletterad med en olineär dämpning som begränsar amplituden) användas som en oscillator.
Lösning: Vi kan räkna med normerade frekvenser, RC
= I.
Då blir
s3 H(s)------ s 3 + 6s 2 + 5s + 1
och problemet är att avgöra när det finns lösningar i höger halvplan till ekvationen H(s)
1
= -k.
Enligt argumentvariationsprincipen är antalet lösningar lika med antalet poler plus antalet varv kring -1/k, så vi behöver bestämma antalet poler och Nyqvistkurvan. Av fysikaliska skäl är det klart att det öppna systemet är stabilt (passivt RCnät), och alltså att antalet poler till H i höger halvplan är noll. Vi skall snart göra sådana undersökningar matematiskt men övergår nu till att rita Nyqvistkurvan. Denna är bilden av positiva imaginära axelns= iw. w > 0. Men . H ( iw)
= -
(iw) 3 -iw 3 ("iw )3 + 6("iw )2 + oiw ~ . + I = -iw . 3 - &,,2 + oiw ~. +I -iw 3 ((1 - 6w 2 ) - iw(5 - w 2 )) (1 - 6w 2 ) 2 + (w(5 - w 2 ))2 ·
KAPITEL 15. STABILITETSTEORI
336
Alltså är
där N blir w
. ) Re H( iw
= w4 (w 2 -
. ) I m H( zw
= wl(&,,2 -1)
N
5)
.V
= (1-6w2)2+ (w(5-w 2)) 2. De kritiska frekvenserna. som ger axelskärningar, = 0, w = 1/ y'6 och w = Js. Vi får följande värdetabell för Nyqvistkurvan: 0 ReH(iw) 0 Im H(iw) 0 w
1/v'6 v'5 0 -1/29 0 + i 5v'5/29
R
+ 1 + 0
Med hjälp av tabellen kan vi skissera Nyqvistkurvan och H(C). Till vänster ser vi en grov skiss gjord för hand efter tabellen ovan, medan den högra bilden visar en mer exakt beräkning med Matlab. Observera här en skillnad mellan konventionerna i matematiken och i reglertekniken. I matematiska sammanhang räknas kurvor som positivt orienterade om de går runt motsols (moturs regniga dagar). Detta medför att i vår standardsituation med halvcirkeln såg~ nomlöps Nyqvistkurvan i riktning från höga frekvenser till låga. I reglertekniken sätter man däremot ut pilar i riktning av växande w, och orienteringen blir som i den högra kurvan.
1 29
J
~~. i :',,
1, :
J
:.
~
-----,E-----
.
0
--
UQAUU
Lägg märke till att det är minst lika lätt att räkna varv i den grova skissen (och det krävs mycket goda ögon för att se att den dubbla axelskärningen ligger i -1/29). Varvräkning i någon av figurerna ger, eftersom antalet poler P till det öppna syste-
15.6. ROUClmS SATS OCH ALGEBRANS FUNDAMENTALSATS
337
met är noll: antal varv kring w 0 a >1 1 O 6. Föra = 2 finns ett nollställe (z = 0) på imaginära axeln, föra = 6 två stycken (z = ±2i). D Metoden ovan kan också användas för att undersöka stabilitet av tillståndsekvationer. Systemet dx dt =Ax+ f(t)
15.7. POLER PÅ RANDEN OCH I OÅNDLIGHETEN
341
är ju stabilt precis då alla egenvärden till A, alltså alla lösningar till polynomekvationen PA(s) = det(s/ - A) = 0 ligger i vänster halvplan. Vi kan alltså undersöka stabilitet av A genom att använda argumentprincipen på dess karakteristiska polynom. För system av låg ordning kan denna metod ibland vara lämplig vid handräkning, även om i sådana fall Routh-Hurwitz kriterier brukar vara bekvämare att använda. För system av högre ordning är det av numeriska stabilitetsskäl i allmänhet olämpligt att beräkna det karakteristiska polynomet, och den metod man får använda är direkt numerisk bestämning av egenvärdena till A. Vi övergår nu till ett annat fall där också obegränsade funktioner förekommer. nämligen då den undersökta funktionen F(s) har en pol på randkurvan C. Då kan argumentprincipen inte användas direkt, ty den förutsätter ju att F( s) är analytisk på C. Det är dock lätt att modifiera metoden så att den fungerar även i dessa fall. Antag t ex att vi är intresserade av höger halvplan och att F har en pol av ordning v i punkten s = p = iw0 • Vi kringgår då singulariteten med en liten halvcirkel med centrum i p och radie r. och låter alltså C vara en i kanten något naggad stor halvcirkel. Vilket bidrag till argumentvariationen ger då en liten halvcirkel kring en pol till argumentvariaP tionen och till bildkurvan C? Eftersom polens ordning är v. så kan vi skriva F(s) = (s - p)- 11 G(s) där G(s) är analytisk is= p och G(p) -I 0. Den lilla cirkeln kan parametriseras 7r
7r
--