Kapitel 26
Integraltransformationen und Z-Transformation
Maple stellt im Paket inttrans folgende Transformationen zur Verfügung:
> with(inttrans);
Das sind also drei Fouriertransformationen, sowie die Laplace-, Mellin-, Hankel- und Hilberttransformation. Mit 'addtable' und 'savetable' können zu den jeweiligen Transformationen Namensgebungen und Rechenregeln hinzugefügt und abgespeichert werden. Mit ?examples,index finden Sie Beispiele unter 'Integral Transforms'.
Die Fouriertransformation wurde wegen ihrer herausragenden Bedeutung (und als Bindeglied Reihenentwicklung - Integraltransformation) schon im vorangehenden Kapitel behandelt.
In diesem Kapitel steht bei den Integraltransformationen die Laplacetransformation im Vordergrund. Im Zeitalter der Digitalisierung darf aber ihre diskretisierte Verwandte, die Z-Transformation, nicht fehlen. (Die Z-Transformation befindet sich nicht im Paket inttrans.)
>
Bei der Laplacetransformation (kurz L-Transformation) wird vorausgesetzt, dass die Originalfunktion f(t) für t<0 Null ist und für t>=0 stückweise glatt. Die Laplacetransformation ist folgendermaßen definiert:
> restart:with(inttrans):
> F(s):=convert(laplace(f(t), t, s),int);
Im Gegensatz zur Fouriertransformation ist der Exponent des Kerns nun nicht mehr rein imaginär, sondern s ist eine komplexe Zahl (mit einem Realteil), es kann also Probleme mit der Konvergenz des Integrals geben. Wenn die Originalfunktion wie exp(a*t) wächst, so ist die Bildfunktion nur für Re(s)>a definiert. Maple kann dem Benutzer die nun erforderlichen Fallunterscheidungen nicht abnehmen. Ein einfaches Beispiel mag dies verdeutlichen:
> f:=t->1;
> value(F(s));
Definite integration: Can't determine if the integral is convergent.
Need to know the sign of --> s
Will now try indefinite integration and then take limits.
> assume(s>0):%;
Bitte behalten Sie dies im Hinterkopf, wenn Maple entweder kein Ergebnis liefert, oder nicht das erwartete. Man hört immer wieder (oder liest es in den Mailinglisten): 'Warum vereinfacht Maple den Term nicht so wie ich mir das vorstelle?'. In den meisten Fällen kann das aber darauf zurückgeführt werden, daß der Benutzer in reellen Zahlen denkt oder sich um die Konvergenz eines Integrals nicht kümmert, während Maple für solche Fälle Stoppstellen eingebaut hat, über die man nur kommt, wenn man die passenden Annahmen macht (der Dialog zwischen CAS und Benutzer läßt sich mit ?infolevel fördern).
Nun ein paar Beispiele zum Gebrauch von 'laplace' und 'invlaplace':
> s:='s':
> laplace(1,t,s);
> invlaplace(%,s,t);
> laplace(t^n,t,s);
Definite integration: Can't determine if the integral is convergent.
Need to know the sign of --> n+1
Will now try indefinite integration and then take limits.
So ähnlich hatten wir das doch schon einmal?
> assume(n,posint):
> p:=laplace(t^n,t,s);
> invlaplace(p,s,t);
Testen Sie auch andere Annahmen wie z.B. assume(n,negint): und assume(n,nonposint):
Wegen der Voraussetzung, dass die Originalfunktion f(t) für t<0 verschwindet, spielt die Heavisidefunktion im Zusammenhang mit der L-Transformation eine wichtige Rolle
>
> alias(sigma=Heaviside):
> laplace(sigma(t),t,s);
Das ist die Transformierte von 1, weil Maple davon ausgeht, dass der Benutzer mit 1 meint: 'f(t)=1 für t>=0 und Null für t<0'. Die verschobene Sprungfunktion erfordert aber bereits eine Fallunterscheidung:
> a:='a':laplace(sigma(t-a),t,s);
> assume(a>0):%;
> invlaplace(%,s,t);
> assume(a<0):%%%;
> invlaplace(%,s,t);
Die Verschiebung nach links geht also bei der Rücktransformation verloren, weil die Originalfunktion laut Annahme immer links von Null verschwindet.
>
Weitere Beispiele zur Verschiebung von Funktionen
So wird eine Funktion korrekt verschoben (Verschiebungssatz):
> f:='f':versch:=laplace(f(t-a)*sigma(t-a),t,s);
> assume(a>0):%;
> assume(a<0):%%;a:='a':
>
So wird eine Funktion abgeschnitten:
> abschn:=laplace(f(t)*sigma(t-a),t,s);
> assume(a>0):%;
> assume(a<0):%%;a:='a':
>
So wird meistens nicht das erreicht, was man will:
> linears:=laplace(f(t-a),t,s);
> assume(a>0):%;
> assume(a<0):%%;a:='a':
>
Konkretes Beispiel
> f:=x->x^2; # hier können eigene Funktionen getestet werden
> versch;
D.h., es wird zunächst nach dem Linearitätssatz ausgewertet
> assume(a>0): simplify(versch);
In Übereinstimmung mit dem Verschiebungssatz (Verschiebung nach rechts)
> assume(a<0):simplify(versch);
In Übereinstimmung mit dem Verschiebungssatz (Verschiebung nach links, oben ist a<0)
> simplify(exp(a*s)*(laplace(f(t),t,s)-int(exp(-s*t)*f(t),t=0..a)));
>
Abschneiden der unverschobenen Funktion rechts vom Ursprung ist keine Verschiebung nach links:
> assume(a>0): simplify(abschn);
und hat links vom Ursprung wegen der Annahme f(t)=0 für t<0 keine Wirkung
> assume(a<0):simplify(abschn);
>
Eine Verschiebung der Originalfunktion ohne Berücksichtigung von f(t)=0 für t<0 ist nicht im Sinne des Verschiebungssatzes
> assume(a>0): simplify(linears);
was aber bei einer Verschiebung nach links keinen Fehler verursacht
> assume(a<0):simplify(linears);
>
Die zwei wichtigsten Eigenschaften der L-Transformation
1. Die Transformierte von Ableitungen (und Integralen)
> restart:with(inttrans):
> alias(F(s)=laplace(f(t),t,s)):
> laplace( diff(f(t), t), t, s);
> laplace( diff(f(t), t$2), t, s);
> laplace(int(f(T), T=0..t), t, s);
D.h., die L-Transformierte der Ableitung einer Funktion ist s mal die L-Transformierte der Funktion (minus Anfangswert). Damit können Differentialgleichungen algebraisiert werden:
>
> dgl:=diff(f(t), t$2)+ diff(f(t), t)+f(t)=0;
Transformation
> alggl:=laplace( dgl ,t,s);
Lösung im Bildbereich
> bildsol:=solve(alggl,F(s));
Rücktransformation
> sol:=invlaplace(bildsol,s,t);
Das läuft also ab, wenn man dsolve mit method=laplace aufruft:
> dsolve(dgl,f(t),method=laplace);
Zum Vergleich
> dsolve(dgl,f(t));
Vordergründig besteht der Vorteil der Laplacemethode gegenüber 'dsolve-pur' also darin, dass man die Anfangsbedingungen 'frei Haus' bekommt und sich nicht mit den Integrationskonstanten herumschlagen muß. Der tiefere Sinn der Laplacemethode liegt aber in der Anpassung der Methode an das Problem. Diese Anpassung kann noch weiter verfeinert werden, wenn der Benutzer die Methode nicht über dsolve aufruft, sondern obige Schritte selbst durchführt, um im Bedarfsfall Maple helfen zu können.
>
2. Faltungssatz
> restart:with(inttrans):
Zur besseren Übersicht werden folgende Abkürzungen verwendet:
> alias(L=laplace):alias(F1(s)=laplace(f1(t),t,s)):alias(F2(s)=laplace(f2(t),t,s)):
>
Zur L-Transformation gehört die einseitige Faltung (von 0 bis t)
> faltung:=int(f1(t-tau)*f2(tau),tau=0..t);
> bild:=L(faltung,t,s);
In Worten: Die Transformierte einer Faltung ist das Produkt der Transformierten - Algebraisierung!
>
Auch die Umkehrung ist Maple bekannt:
> invlaplace(bild,s,t);
>
Anwendung siehe Lösung von DGLn mit der L-Transformation
>
Für Faltungen wie sie bei Messungen vorkommen ist die L-Transformation weniger geeignet, weil die L-Faltung einseitig ist. Beispiele zur zweiseitigen Faltung finden Sie bei der Fouriertransformation.
>
Addtable, savetable
Zu addtable finden Sie einige Beispiele in der Maple-Hilfe: ?examples/addtable
Hier ist noch ein Beispiel wie man neben Namensgebungen und Rechenregeln konkrete Einträge hinzufügt:
> restart:with(inttrans):
> laplace(exp(-x^2),x,s);
> invlaplace(%,s,x);
Die Rücktransformierte ist zunächst nicht bekannt, kann aber mit addtable hinzugefügt werden
> addtable(invlaplace,exp(1/4*s^2)*erfc(1/2*s),2/sqrt(Pi)*exp(-x^2),s,x);
Mit infolevel kann man Maple beim Suchen zuschauen
> infolevel[all]:=5:
> invlaplace(1/2*sqrt(Pi)*exp(1/4*s^2)*erfc(1/2*s),s,x);
inttrans/invlaplace/lookup: Entered table lookup (invlaplace) exp(1/4*s^2)*erfc(1/2*s)
inttrans/invlaplace/lookup: Looking up exp(1/4*s^2)*erfc(1/2*s) in Whittaker table
inttrans/invlaplace/lookup: Looking up exp(1/4*s^2)*erfc(1/2*s) in other table
inttrans/invlaplace/lookup: Returning from lookup (invlaplace) 2/Pi^(1/2)*exp(-x^2)
Beim zweiten Aufruf wird die Information nicht mehr angezeigt
> invlaplace(1/2*sqrt(Pi)*exp(1/4*s^2)*erfc(1/2*s),s,x);
weil der Eintrag nun auch in der remember table von invlaplace steht (dort schaut Maple zuerst nach)
> op(4,eval(invlaplace));
>
Bitte beachten Sie, daß der Eintrag ohne den konstanten Vorfaktor gemacht werden muß:
> restart:with(inttrans):
> laplace(exp(-x^2),x,s);
> invlaplace(%,s,x);
> addtable(invlaplace,1/2*sqrt(Pi)*exp(1/4*s^2)*erfc(1/2*s),exp(-x^2),s,x);
>
> infolevel[all]:=2:
> invlaplace(1/2*sqrt(Pi)*exp(1/4*s^2)*erfc(1/2*s),s,x);infolevel[all]:='infolevel[all]':
inttrans/invlaplace/lookup: Entered table lookup (invlaplace) exp(1/4*s^2)*erfc(1/2*s)
simplify: applying commonpow function to expression
simplify: applying power function to expression
combine: combining with respect to combine/exp
combine: combining with respect to combine/exp
simplify: applying simplify/power function to expression
inttrans/invlaplace/lookup: Entered table lookup (invlaplace) exp(1/4*tt^2)
combine: combining with respect to combine/exp
simplify: applying simplify/power function to expression
inttrans/invlaplace/lookup: Entered table lookup (invlaplace) erfc(1/2*tt)
simplify: applying simplify/power function to expression
inttrans/invlaplace/lookup: Entered table lookup (invlaplace) exp(1/4*tt^2)*erfc(1/2*tt)
>
Die Tabelleneinträge können mit savetable gespeichert und wieder eingelesen werden werden.
>
>
Lösung von Differentialgleichungen mit der Laplacetransformation
Sprungantwort eines RC-Glieds
Das prinzipielle Vorgehen beim Lösen von DGLn wurde schon beschrieben. Hier sind weitere konkrete Beispiele.
Die Sprungantwort eines RC- oder RL-Glieds kann durch eine DGL erster Ordnung beschrieben werden
> restart:with(inttrans): alias(sigma=Heaviside):
> de:=diff(y(t),t)+b*y(t)=u;
> u:=sigma(t-1)-sigma(t-5);
Transformsation der DGL
> alias( Y(s)=laplace(y(t),t,s) ):
> laplace(de, t, s);
Auflösen nach der Transformierten
> sol:=solve(subs(y(0)=0,%), Y(s));
Rücktransformation
> antwort:=invlaplace(sol ,s,t);
>
> b:=2:
> plot({b*antwort,u}, t=0..9);
>
RLC-Kreis
Ebenfalls eine bekannte DGL, diesmal explizit mit der L-Transformation gelöst
> restart:with(inttrans):
> sys:= u(t)=R*j(t) + L*diff(j(t),t) + 1/C*int(j(t), t):
> de:= diff(sys, t);
> alias( U(s)=laplace(u(t), t, s), J(s)=laplace(j(t), t, s)):
Transformation
> bild:=laplace( de, t, s);
Auflösen nach dem transformierten Strom
> Js:=solve(bild, J(s));
Einsetzen von Anfangsbedingungen
> Js:=subs( u(0)=0, j(0)=0, D(j)(0)=0, Js);
Bereitstellen einer transformierten Spannung (hier im Originalbereich konstant)
> U(s):=laplace(ut, t, s);
> jt:=invlaplace(Js, s, t);
Werte einsetzen
> R:=12: L:=2: C:=0.02: ut:=200:
> plot({ut, jt*R}, t=0..2);
Sinusförmige Erregerspannung
> R:='R':C:='C':L:='L':U(s):='U(s)':
> ut:=100*sin(2*Pi*50*t):
> us:=laplace(ut, t, s);
> invlaplace( subs( U(s)=us, Js), s, t):
> jt:=subs(R=10, L=2, C=0.02, %): simplify(%): evalf(%,5);
> plot( {jt*10}, t=0..1, numpoints=2000);
Trapezförmige Erregerspannung
> alias(sigma=Heaviside):
> ut:=t-(t-1)*sigma(t-1)-(t-4)*sigma(t-4)+(t-5)*sigma(t-5);
> us:=laplace(ut, t, s);
Werte einsetzen
> subs( R=8, L=2, C=0.02,U(s)=us, Js);
> invlaplace(%,s,t): jt:=evalf(%,5);
> plot( {ut, jt*16}, t=0..8);
>
Wärmeleitung (partielle DGL)
Eindimensionale Wärmeleitungsgleichung (homogenes Medium, b ist der Leitfähigkeit proportional):
> restart:with(inttrans):
> pde:=b^2*diff(u(t,x),x$2)-diff(u(t,x),t)=0;
Lösungsversuch mit pdsolve
> sol:=pdsolve(pde,u);
> solb:=pdsolve(pde,u,build);
Kann physikalisch nicht sein: _c[1] müsste negativ sein, das würde aber eine Oszillation mit x bedeuten.
L-Transformation
> alias(U=laplace(u(t,x),t,s)):
> bildgl:=laplace(pde,t,s);
Für die Lösung der DGL benötigt dsolve eine Funktion von x, die Anfangsbedingung wird 0 gesetzt (kann nach der Rücktransformation wieder hergestellt werden).
> de:=subs(U=T(x),u(0,x)=0,bildgl);
Transformierte Randbedingungen A(s)
> bildsol:=dsolve({de,T(0)=A0(s),T(l)=A1(s)},T(x));
Ranbedingung mit l -> infinity
> bildsoll:=limit(bildsol,l=infinity);
Kann nur ausgewertet werden, wenn über s und b nähere Angaben gemacht werden
> assume(s>0,b>0):bildsolla:=bildsoll;
Maple benötigt für den Faltungssatz noch die Angabe, dass A0(s) eine transformierte Randbedingung ist
> A0(s):=laplace(a0(t),t,s):
> sol:=invlaplace(rhs(bildsoll),s,t);
Angabe einer Randbedingung im Originalbereich
> subs(a0(_U1)=a0,sol);
> slov:=value(%);
Auch hier wieder eine nähere Angabe
> assume(x>0):solv:=%;
Nun müssen die Annahmen wieder entfernt werden, damit mit der Lösung weitergearbeitet werden kann, diesmal mit einem Trick:
> solv:=parse(convert(solv,string));
Konkrete Werte
> b:=0.5:a0:=1:
> plot3d(solv,t=0.001..10,x=-5..5,axes=framed,grid=[10,50],style=patchcontour,contours=20,shading=Z);
>
Weitere Integraltransformationen
Hier werden zur Vollständigkeit alle Transformationen aufgelistet, die in Maple bekannt sind (außer der Laplacetransformation), zusätzlich wird jeweils die Rücktransformation und die Transformierte von angegeben (sofern verfügbar). Beispiele finden Sie in ?examples,index -> Integral Transforms.
Die Fouriertransformation
> restart: with(inttrans):
> convert(fourier(f(t),t,s),int);
und ihre Rücktransformation
> convert(invfourier(f(t),t,s),int);
> fourier(diff(f(t),t),t,s);
>
Fouriersinus
> convert(fouriersin(f(t),t,s),int);
ist selbstinvertierend
> fouriersin(fouriersin(f(t),t,s),s,t);
> fouriersin(diff(f(t),t),t,s);
>
Fouriercosinus
> convert(fouriercos(f(t),t,s),int);
> fouriercos(fouriercos(f(t),t,s),s,t);
> fouriercos(diff(f(t),t),t,s);
(enthält eine Anfangsbedingung)
>
Mellintransformation
> convert(mellin(f(t),t,s),int);
Zur Rücktransformation siehe ?examples,mellin
insbesondere:
"The lesson to be learned from the above is that you must have the correct assumptions on parameters for all integral transforms, and the correct range for the inverse Mellin transform."
>
Hankeltransformation
> convert(hankel(f(t),t,s,v),int);
> v:='v':
> hankel(hankel(f(t),t,s,v),s,z,v);
> assume(v>-1/2,v<>0):%;
>
> hankel(diff(f(t),t),t,s,v);
>
Hilberttransformation
> convert(hilbert(f(t),t,s),int);
> invhilbert(f(s),s,t);
> hilbert(hilbert(f(t),t,s),s,t);
> hilbert(diff(f(t),t),t,s);
>
Die Z-Transformation steht im Maple-Browser zwar unter Integraltransformationen, ist aber eine diskrete Transformation und kann auch angesprochen werden, ohne das Paket inttrans zu laden.
Einführung
Die Z-Transformation
> restart:
> F(z) = sum(f(n)/z^n,n=0..infinity);
entsteht aus der L-Transformation durch Diskretisierung (Digitalisierung):
Tastet man stetige Funktionen f(t) mit einer Samplingrate a ab, so entstehen z.B. Treppenfunktionen
> digi:=(f,a)-> f(floor(a*t)/a);
> plot([10*sin(t),10*digi(sin,5),digi(x->x^2,10)/5],t=-2..7,numpoints=500);
> #plot(sin(x),x=0..10,adaptive=false,sample=[seq(i,i=0..10)]); #so nicht :-))
>
Auf diese Treppenfunktionen läßt sich nun die L-Transformation
> F(s)=convert(laplace(f(t),t,s),int);
folgendermaßen anwenden: die Zeit t wird diskretisiert zu n*T, die L-Transformation innerhalb T für f(n)=const durchgeführt und aufsummiert, also z.B. für T=1
>
> sum(Int(f(n)*exp(-s*t),t=n..n+1),n=0..infinity)=sum(int(f(n)*exp(-s*t),t=n..n+1),n=0..infinity);
oder (mit der Maus :-) vereinfacht:
> L:=(1-exp(-s))/s*sum(f(n)*exp(-s*n),n = 0 .. infinity);
Dabei nennt man
> DL:=sum(f(n)*exp(-s*n),n = 0 .. infinity);
die diskrete L-Transformation. Setzt man hierin exp(s)=z, so wird daraus
> subs(exp(-s*n)=1/z^n,DL);
und zwischen der Z-Transformation und der L-Transformation besteht die Korrespondenz
> subs(exp(-s)=1/z,L/DL)*subs(exp(-s*n)=1/z^n,DL)*s=s*convert(laplace(f(t),t,s),int);
oder
> (1-1/z)*F(z)=s*F(s);
Korespondenz heißt: es muss eine Treppenfunktion sein:
> (1-1/z)*ztrans(n,n,z)=s*inttrans[laplace](floor(t),t,s);
> simplify(%);
> simplify(subs(exp(-s)=1/z,%));
und nicht eine stetige Funktion:
> simplify((1-1/z)*ztrans(n,n,z)=s*inttrans[laplace](t,t,s));
Bei der Verwendung von ztrans ist zu beachten, dass z und F(z) im allgemeinen Fall komplex sind, und dass die Transformation nur definiert ist, wenn die Reihe konvergiert. Einfachstes Beispiel:
> ztrans(a,n,z);
Das ist der Grenzwert der geometrischen Reihe
> Sum(a/z^n,n=0..infinity)=sum(a/z^n,n=0..infinity);
der nur für |1/z|<1 existiert.
>
Die Rücktransformation gibt es natürlich auch:
> invztrans(rhs(%),z,n);
Weitere Beispiele
> ztrans( (-1)^n, n, z);
> simplify(ztrans(a^n,n,z));
> invztrans(%,z,n);
Einen Sonderfall (bzw. eine Analogie zur Diracfunktion) stellt die 'charakteristische Funktion' dar, die man erhält, wenn man eine Konstante im Bildbereich zurücktransformiert:
> invztrans(a,z,n);
Dabei meint eine Folge, die nur für n=0 1 ist und Null sonst:
> seq(charfcn[0](n),n=-3..3);
und es gilt z.B.
> invztrans(1/z^5,z,n);
> ztrans(f(n)*charfcn[5](n),n,z);
d.h., es wird aus der unendlichen Reihe der 5. Summand isoliert.
(Die charakteristische Funktion kann auch zum Morsen verwendet werden:
> seq(charfcn[-1,1,2](n),n=-3..4);
> seq(charfcn[n](n),n=-3..4);
> seq(charfcn[n mod 2](2*n mod 4),n=-3..20);
aber das nur nebenbei)
>
In Analogie zur L-Transformation besitzt die Z-Transformation die Eigenschaft, dass die Transformierte der Folge der Differenzen im wesentlichen das (z-1)-fache der Transformierten der Folge ist:
Erste Differenzen
> d1:=n->f(n+1)-f(n);
> alias(F(z)=ztrans(f(n),n,z)):
> collect(ztrans(d1(n),n,z),F(z));
Auf dieser Eigenschaft bauen die Lösungsverfahren für Differenzengleichungen auf
Kleiner Exkurs zu höheren Differenzen
Zweite Differenz
> d2:=n->d1(n+1)-d1(n);
> d2(k);
> collect(ztrans(d2(n),n,z),F(z));
Man kann das automatisieren:
a)
> d[1]:=n->f(n+1)-f(n):
>
for i from 2 to 5
do d[i]:=unapply(d[i-1](n+1)-d[i-1](n),n); od:
> d[5](k);
b)
> dn:=n->add((-1)^(i+n)*binomial(n,i)*f(k+i),i=0..n);
> dn(5);
c)
>
differ:=proc(n,k)
if n<2 then f(k+n)-f(k+n-1);
else differ(n-1,k+1)-differ(n-1,k) end if
end proc;
>
> differ(5,k);
> differ(5,1);
> seq(differ(i,k),i=1..5);
Wenn z.B. die zweite Differenz konstant ist,
> ztrans(differ(2,n)=a,n,z);
>
> solve(%,F(z));
> collect(invztrans(%,z,n),n);
so muss die Folge quadratisch wachsen.
Geht man von den Differenzen zu den Differenzenquotienten über (Division durch die Schrittweite h^n), so kann man Differentialgleichungen durch Differenzengleichungen annähern, wichtigste Anwendung: finite Elemente.
>
>
Lösung von Differenzengleichungen mit der Z-Transformation
Erstes Beispiel
Die 'Schwingungsgleichung'
> alias(F(z)=ztrans(f(k),k,z)):
> dgl:=f(k+2)+f(k)=0;
wird Z-transformiert
> bild:=ztrans(dgl,k,z);
> bildsol:=solve(bild,F(z));
nach der Transformierten aufgelöst und zurücktransformiert
> sol:=invztrans(bildsol,z,k);
Lassen Sie sich nicht abschrecken - das funktioniert!
> sol:=eval(sol,{f(0)=0,f(1)=5});
Man kann es sogar zeichnen:
> plot(sol,k);
Was haben wir da gezeichnet?
> add(allvalues(sol)[i],i=1..2);
> evalc(Re(%));
Wir werden das in einem zweiten Beispiel näher untersuchen, zunächst der Vergleich mit rsolve:
> evalc(Re(rsolve(dgl,f)));
>
Zweites Beispiel
Paratstellen der Digitalisierungsfunktion und höherer Differenzen
> restart:
> digi:=(f,a)-> f(floor(a*k)/a):
>
differ:=proc(n,k)
if n<2 then u(k+n)-u(k+n-1);
else differ(n-1,k+1)-differ(n-1,k) end if
end proc:
> alias(Delta=charfcn[0], U(z)=ztrans(u(k),k,z)):
>
Die erste und zweite Differenz stehen nun als Differenzenquotienten mit der Schrittweite 1 für die erste und zweite Ableitung. Es können verschiedene Anregungen ('rechte Seiten') untersucht werden. Der schwarze Input läßt sich mit der rechten Maustaste zu 'Mapleinput' toggeln.
1. Ungestörtes System
> dgl:=8*differ(2,k)+4*differ(1,k)+3*u(k)=0;
2. Einmalige Auslenkung
dgl:=8*differ(2,k)+4*differ(1,k)+3*u(k)=100;
3. Impuls zur Zeit Null
dgl:=8*differ(2,k)+4*differ(1,k)+3*u(k)=100*charfcn[0](k);
4. Sinusförmige Anregung
dgl:=8*differ(2,k)+4*differ(1,k)+3*u(k)=100*sin(2*Pi*k/3);
Transformation der Differenzengleichung und Auflösen nach der Bildfunktion
> bild:=ztrans( dgl, k, z);
> Uz:=solve(bild, U(z));
Rücktransformation
> uk:=invztrans( Uz, z, k);
>
Wer sich für _alpha interessiert, kann hier nachsehen
allvalues( op(selectfun(uk,RootOf)));
Zum Vergleich mit rsolve
rsolve(dgl,u);
>
Anfangswerte
> uka:=subs( u(0)=100, u(1)=50, uk);
Kontrollausgabe
> seq(uka, k=0..10); #k:='k':
Plot
> plot(uka, k=0..20);
Nanu - haben wir nun eine stetige Funktion?
> plot(uka, k=0..10,adaptive=false,sample=[seq(i,i=0..10)]);
> ukf:=unapply(uka,k);
> plot([ukf(k),digi(ukf,1)],k=0..20);
Also doch digital!
>
Neue Anfangswerte, neue Differenzengleichung
>