White-Box-Systemidentifikation automatisieren

[Tobias S. (tobi_s)]

Hallo!

Ich lese oft im Forum hier mit und habe mich nun auch mal angemeldet, 
weil ich eine Frage habe! ;)

Ich bin momentan mit der Evaluation verschiedener 
Systemidentifikationsverfahren für für lineare- und nichtlineare Systeme 
beschäftigt. In der Literatur bin ich dabei auf sehr präzise Algorithmen 
gestoßen, welche aus gemessenen Eingangs- und Ausgangssignalen eine 
Aussage über die Systemstruktur- und Parameter machen können.

Dabei ist es mir gelungen aus Simulationsdaten von Simulatoren (SPICE / 
Simulink Simscape) Blackbox-Identifikationen vorzunehmen. Die Systeme 
führe ich in ein diskretes Zustandsraummodell über, um sie später auf 
digitalen Signalprozessoren implementieren zu können.

Da ich oft mit Systemen zu tun habe, deren innere Parameter variabel 
sind - sprich elektronische Schaltungsmodelle - ist der Black-Box Ansatz 
natürlich nicht immer optimal. Verschiedene Matrizen über 
Lookup-Tabellen zu laden - je mehr variable Bauelemente, umso mehr 
Dimensionen ergeben sich, also exponentielles Wachstum der Daten.

Deshalb habe ich einige Experimente in der symbolischen Analyse von 
elektrischen Netzwerken vorgenommen, wobei ich in der amerikanischen 
Literatur auch auf das nichtlineare Zustandsraummodell gestoßen bin. Der 
Vorteil der symbolischen Formulierung besteht natürlich darin, dass ich 
bei Vorgabe der Bauelemente sehr schnell die notwendigen 
Matrizenparameter berechnen kann (Bei LTI lediglich 4 Matrizen, A, b, c, 
d). Nun ergibt sich für mich folgendes Problem: Die Schaltungsanalyse 
ist oft ein sehr langwieriger Prozess. Zwar ist es häufig möglich eine 
Schaltung in viele Subsysteme zu untergliedern, doch ziehen bereits 
kleine Änderungen der Schaltung große Auswirkungen in der mathematischen 
Formulierung nach sich. Leistungsfähige CAS-Software erleichtert mir die 
Berechnungen zwar immens, doch alle mir bekannten Simulatoren arbeiten 
(sinnvollerweise) rein numerisch.

Nun meine Frage: Gibt es eine Softwarelösung, welche mir aus einer 
gegebenen Netzliste (oder sonstiger Schaltplaneingabe...) einen Satz 
mathematischer Gleichungen gemäß den Kirchhoffschen Gesetzen liefern 
kann? Beispielsweise zur Weiterverarbeitung in MATLAB/Simulink [-> 
Überführung in diskretes Zustandsraummodell]. Momentan führe ich die 
Berechnugnen "von Hand" durch, was extrem zeitraubend ist. Der Prozess 
stellt sich für mich sehr algorithmisch dar, weshalb ich mir nicht 
vorstellen kann, dass keine Softwarelösung existiert? Man bedenke nur 
die Berücksichtigung von nichtidealen Eigenschaften von Bauelementen 
oder umfangreicheren Transistormodellen.

Hier ein einfaches Beispiel einer Modellierung eines reinen LTI-Systems: 
http://www.dafx.ca/proceedings/papers/p_001.pdf . Es ergeben sich eben 
häufig "länglich" Ausdrücke der Koeffizienten und wenn man solche 
Diskretisierungsverfahren auf verschiedene Schaltungen - später 
deutliche komplexere - anwenden möchte, dazu ein großer Aufwand. 
Glücklicherweise kann CAD-Software diese Rechnerei bei gegebenen 
Gleichungen übernehmen.

Auf Seite 124 der Dissertation von David Yeh 
(https://ccrma.stanford.edu/~dtyeh/papers/DavidYehThesissinglesided.pdf) 
ist ein Beispiel dieser nichtlinearen Zustandsraummodellierung einer 
analogen Schaltung gegeben. Diesen Prozess (Herleitung der Matrizen) 
möchte ich letzendlich automatisieren.

Mir würde es zunächst genügen, eine Software(-routine) für die Ausgabe 
der Gleichungen für ein reines LTI-System zu benutzen. Natürlich soll 
die Software nicht die Matrizen herleiten, es genügt mir eine linear 
unabhängige Formulierung der Netzwerkgleichungen. Die Software selbst zu 
schreiben wäre für mich ein weiteres Subprojekt, welches viel Zeit in 
Anspruch nehmen würde, weshalb ich gerne auf eine bereits vorhandene 
Softwarelösung zurückgreifen würde. Für die Weiterverarbeitung 
(Diskretisierung, Lösen von Nichtlinearitäten etc.) würde ich mir eine 
MATLAB-Lösung implementieren.


Ich bedanke mich im Voraus!

[A. S. (rava)]

schau dir mal Dymola an. Das ist ein Simulationstool, vergleichbar mit 
Simulink. Allerdings werden im Design keine Eingänge und Ausgänge 
definiert, man kann "in den Analyser geben" was man möchte und das 
System stellt die Gleichungen intern nach der gesuchten Größe um, falls 
möglich.

Dazu wird ein einheitlicher systemtheoretischer Ansatz verwendet, 
bestehend aus Flussgrößen (z.B. Strom) und Potentialgrößen (z.B. 
Spannung). In der Elektronik ist das noch relativ überschaubar. Die 
zugrundeliegende Modellierungssprache heißt Modelica und da fängst du 
am besten zu suchen an, denn es gibt neben Dymola jede Menge Software, 
die darauf aufbaut.

Das Wolfram-Tool kann z.B. ganz sicher Gleichungssysteme generieren und 
dan Bedarf umstellen. MapleSim sieht auch gut aus. Es kann nämlich sein, 
dass Dymola nur mit numerischen Bauteilwerten arbeitet.
Ist halt auch eine Preisfrage ;)

Ich habe mit keinem der Tools wirklich intensiv gearbeitet und ich weiß 
auch gar nicht, ob es das ist, was du suchst.

hier geht's los:
http://en.wikipedia.org/wiki/Modelica


halt uns auf dem Laufenden!
klick: https://www.youtube.com/watch?v=s8J8nvt7ZAA

[Tobias S. (tobi_s)]

Hallo nochmal!

Vielen dank für die schnelle Antwort! Modelica sieht sehr interessant 
aus! Werde ich mir auf jeden Fall zu Gemüte ziehen.

Ich habe soeben ein Softwarepaket aus der Welt der freien Software 
gefunden,  QSapecNG, es handelt sich um ein Frontend für sapec-ng, 
welches Kommandozeilenbasiert mit Netzlisten arbeitet. QSapecNG bietet 
einen komfortablen Editor an, so kann man sich das tippen der Netzlisten 
sparen ;) !

Erstaunlicherweise ist das Programm extrem einfach zu benutzen. 
Eingangsspannungsquelle definieren, Ausgangsport definieren - 1 Klick 
und das Programm berechnet die Laplace-Übertragungsfunktion. Die 
Software scheint wohl exakt für die LTI-Systemanalyse ausgelegt worden 
zu sein. Ich musste nicht einmal spezielle Optionen festlegen! :o

Das ganze vereinfacht meine Arbeit enorm. Vielleicht kann ich mir die 
ganze Prozedur mit Modelica ja noch einfacher und vor allem besser und 
umfangreicher gestalten. Ich bin auch in der Regelungstechnik unterwegs, 
deshalb kommt es mir sehr gelegen eine solch mächtige 
Modellierungssprache zu nutzen. Systemidentifikation mechanischer 
Systeme hilft mir sehr (!) für digitale Regelungen im Zustandsraum 
weiter. Schließlich ist die Modellbildung besonders wichtig, wenn 
Zustandsgrößen geschätzt werden sollen.

Habe mal die einfache Tonestack-Schaltung aus dem Beispiel 
(http://www.dafx.ca/proceedings/papers/p_001.pdf) in den Simulator 
gegeben (Screenshot, wie ihr seht, komm ich aus dem Audiobereich). Jetzt 
wird es ein klein wenig "unübersichtlich" ;) Hier meine ersten 
Resultate:

Symbolisches Analyseergebnis, direkt aus QSapecNG:

(( (1-m)*R3  C1  C2  C3  R4  V2  m*R3 + C1  C2  C3  R4  V2 * 
l*R2 * m*R3 + (1-m)*R3  C1  C2  C3  R4  V2  t*R1 + (1-m)*R3  C1  
C2  C3  V2 * m*R3 * t*R1 + C1  C2  C3  R4  V2 * l*R2 * t*R1 + C1 * 
C2  C3  V2 * l*R2 * m*R3 * t*R1 + (1-m)*R3 * (1-t)*R1  C1  C2  C3  
V2 * m*R3 + (1-t)*R1  C1  C2  C3  V2 * l*R2 * m*R3 ) * s^3 + ( C1 * 
C3  R4  V2 * m*R3 + (1-m)*R3  C1  C3  R4  V2 + (1-m)*R3  C1  C3 
 V2  m*R3 + C1  C3  R4  V2  l*R2 + C1  C3  V2 * l*R2 * m*R3 + C1 
 C3  R4  V2  t*R1 + C1  C3  V2 * m*R3 * t*R1 + (1-t)*R1  C1  C3 
 V2  m*R3 + C1  C2  R4  V2  m*R3 + (1-m)*R3  C1  C2  R4  V2 + 
C1  C2  R4  V2  l*R2 + C1  C2  R4  V2  t*R1 + C1  C2  V2 * 
m*R3 * t*R1 + (1-m)*R3  C1  C2  V2  t*R1 + C1  C2  V2 * l*R2 * 
t*R1 + (1-t)*R1  C1  C2  V2  m*R3 + (1-m)*R3 * (1-t)*R1  C1  C2 * 
V2 + (1-t)*R1  C1  C2  V2  l*R2 + (1-m)*R3  C2  C3  V2  m*R3 + 
C2  C3  V2 * l*R2 * m*R3 ) * s^2 + ( C1  V2  m*R3 + (1-m)*R3  C1  
V2 + C1  V2  l*R2 + C1  V2  t*R1 + C3  V2  m*R3 + C2  V2  m*R3 + 
(1-m)*R3  C2  V2 + C2  V2  l*R2 ) * s)/(( (1-m)*R3  C1  C2  C3  
R4 * m*R3 + C1  C2  C3  R4  l*R2 * m*R3 + (1-m)*R3  C1  C2  C3  
R4 * t*R1 + (1-m)*R3  C1  C2  C3  m*R3 * t*R1 + C1  C2  C3  R4  
l*R2 * t*R1 + C1  C2  C3 * l*R2 * m*R3 * t*R1 + (1-m)*R3 * (1-t)*R1 * 
C1  C2  C3 * R4 + (1-m)*R3 * (1-t)*R1  C1  C2  C3  m*R3 + (1-t)*R1 
 C1  C2  C3  R4 * l*R2 + (1-t)*R1  C1  C2  C3  l*R2 * m*R3 ) * 
s^3 + ( C1  C3  R4 * m*R3 + (1-m)*R3  C1  C3 * R4 + (1-m)*R3  C1  
C3 * m*R3 + C1  C3  R4 * l*R2 + C1  C3  l*R2 * m*R3 + C1  C3  R4 * 
t*R1 + C1  C3  m*R3 * t*R1 + (1-t)*R1  C1  C3 * R4 + (1-t)*R1  C1  
C3 * m*R3 + C1  C2  R4 * m*R3 + (1-m)*R3  C1  C2 * R4 + C1  C2  R4 
* l*R2 + C1  C2  R4 * t*R1 + C1  C2  m*R3 * t*R1 + (1-m)*R3  C1  
C2 * t*R1 + C1  C2  l*R2 * t*R1 + (1-t)*R1  C1  C2 * R4 + (1-t)*R1 * 
C1  C2  m*R3 + (1-m)*R3 * (1-t)*R1  C1  C2 + (1-t)*R1  C1  C2 * 
l*R2 + (1-m)*R3  C2  C3 * R4 + (1-m)*R3  C2  C3 * m*R3 + C2  C3  
R4 * l*R2 + C2  C3  l*R2 * m*R3 ) * s^2 + ( C1 * m*R3 + (1-m)*R3 * C1 
+ C1 * l*R2 + C1 * t*R1 + (1-t)*R1 * C1 + C3 * R4 + C3 * m*R3 + C2 * R4 
+ C2 * m*R3 + (1-m)*R3 * C2 + C2 * l*R2 ) * s + 1)

Vereinfachung mit Matlab Symbolic Toolbox, sortiert nach Potenzen von s:

((C1*C2*C3*R1*R3^2*m + C1*C2*C3*R3^2*R4*m - C1*C2*C3*R1*R3^2*m^2 - 
C1*C2*C3*R3^2*R4*m^2 + C1*C2*C3*R1*R3*R4*t + C1*C2*C3*R1*R2*R3*l*m + 
C1*C2*C3*R2*R3*R4*l*m + C1*C2*C3*R1*R2*R4*l*t - 
C1*C2*C3*R1*R3*R4*m*t)*s^3 + (C1*C3*R3^2*m + C2*C3*R3^2*m - 
C1*C3*R3^2*m^2 - C2*C3*R3^2*m^2 + C1*C2*R1*R3 + C1*C2*R3*R4 + 
C1*C3*R3*R4 + C1*C2*R1*R2*l + C1*C2*R2*R4*l + C1*C3*R2*R4*l + 
C1*C3*R1*R3*m + C1*C2*R1*R4*t + C1*C3*R1*R4*t + C1*C3*R2*R3*l*m + 
C2*C3*R2*R3*l*m)*s^2 + (C1*R3 + C2*R3 + C1*R2*l + C2*R2*l + C3*R3*m + 
C1*R1*t)*s)/((C1*C2*C3*R1*R3*R4 + C1*C2*C3*R1*R3^2*m + 
C1*C2*C3*R3^2*R4*m - C1*C2*C3*R1*R3^2*m^2 - C1*C2*C3*R3^2*R4*m^2 + 
C1*C2*C3*R1*R2*R4*l - C1*C2*C3*R1*R3*R4*m + C1*C2*C3*R1*R2*R3*l*m + 
C1*C2*C3*R2*R3*R4*l*m)*s^3 + (C1*C3*R3^2*m + C2*C3*R3^2*m - 
C1*C3*R3^2*m^2 - C2*C3*R3^2*m^2 + C1*C2*R1*R3 + C1*C2*R1*R4 + 
C1*C3*R1*R4 + C1*C2*R3*R4 + C1*C3*R3*R4 + C2*C3*R3*R4 + C1*C2*R1*R2*l + 
C1*C2*R2*R4*l + C1*C3*R2*R4*l + C2*C3*R2*R4*l + C1*C3*R1*R3*m - 
C2*C3*R3*R4*m + C1*C3*R2*R3*l*m + C2*C3*R2*R3*l*m)*s^2 + (C1*R1 + C1*R3 
+ C2*R3 + C2*R4 + C3*R4 + C1*R2*l + C2*R2*l + C3*R3*m)*s + 1)

-------------------------

Die Richtigkeit der Analyse scheint tatsächlich zuzutreffen!! :) Das 
Teil arbeitet wirklich rein symbolisch, auch die "simulierten" Potis mit 
x*R und (1-x)*r für 0 < x <= 1 werden korrekt berücksichtigt.

Schauen wir uns "b3" an, also den Koeffizienten der 3. Potenz von "s" im 
Zähler:

(C1*C2*C3*R1*R3^2*m + C1*C2*C3*R3^2*R4*m - C1*C2*C3*R1*R3^2*m^2 - 
C1*C2*C3*R3^2*R4*m^2 + C1*C2*C3*R1*R3*R4*t + C1*C2*C3*R1*R2*R3*l*m + 
C1*C2*C3*R2*R3*R4*l*m + C1*C2*C3*R1*R2*R4*l*t - 
C1*C2*C3*R1*R3*R4*m*t)*s^3

Fassen wir noch die "m,l,t-Terme" zusammen und sortieren etwas, so 
ergibt sich:

$$b_3 = l m (C_1 C_2 C_3 R_1 R_2 R_3 + C_1 C_2 C_3 R_2 R_3 R_4) - m^2(C_1 C_2 C_3 R_1 R_3^2 + C_1 C_2 C_3 R_3^2 R_4) + m (C_1 C_2 C_3 R_1 R_3^2 + C_1 C_2 C_3 R_3^2 R_4) + t C_1 C_2 C_3 R_1 R_3 R_4 - t m C_1 C_2 C_3 R_1 R_3 R_4 + t l C_1 C_2 C_3 R_1 R_2 R_4$$

Die automatischen Berechnungen stimmen mit den Resultaten aus dem Paper 
überein. Jetzt ab damit in den diskreten Zustandsraum! ;)