LS oder RLS Algorithmus

[Chris (Gast)]

Guten Abend!

hat schon einmal jemand von euch den LS (least square) oder rekursiven 
least square Algorithmus (Methode der kleinsten Quadrate) in C/C++ oder 
C# implementiet? Oder kennt einer eine Dll die diesen bereit stellt?

Meine Aufgabe ist es eine Regelstrecke zu identifizieren. Ich habe eine 
ASCII Datei mit x-y-Werten. Diese habe ich eingelesen usw. doch den 
Algo. bekomme ich nicht ans laufen.

Hat jemand eine Ahnung wie man sowas macht?

Vielen Dank für eure Hilfe

[Gast (Gast)]

Da [1] isses recht schön beschrieben.

HTH

[1] http://en.wikipedia.org/wiki/Least_squares

[avr (Gast)]

Hallo Chris,

hier bei Arndt-Bruenner gibt es solch eine Berechnung (Java-Script).

http://www.arndt-bruenner.de/mathe/scripts/regr.htm

Deine Wertepaare kannst du als Text mit Copy-Past eintragen.
Werteanzahl nicht vergessen und deine "Wunschformel" wählen.

Hoffe das past.

avr

[volly (Gast)]

Guten Abend Chris,

ich habe eine solche Identifikation mit einer PT2-Strecke ausprobiert. 
Allerdings bin ich mit LS nicht weiter gekommen. Denn die LS bringt Dir 
lediglich die Koeffizienten eines Polynoms raus, damit diese verglichen 
werden können, müssten alle Regelstreckenarten als Polynom vorliegen - 
Keine Ahnung, ob dies überhaupt funktioniert.

Eine recht gute Lösung habe ich mit dem Levenberg Marquardt Algorithmus 
erreicht, dabei wird nämlich die zu vergleichende Funktion bereits im 
Vorfeld definiert:

Beispiel: http://www.messen-und-deuten.de/lmfit/lmfit.html

Aber grundsätzlich ist der Aufwand recht hoch, wenn Du eine unbekannte 
Regelstrecke hast. Du musst alle möglichen Funktionen definieren, die 
Koeffizienten bestimmen und dann die mit dem kleinsten Fehler auswählen.

Hört sich einfacher an, als es ist - dennoch eine schöne Aufgabe.

Gruß
Volly

[volly (Gast)]

Hallo nochmal,

ich habe gerade mal im Handbuch der Regelungstechnik nachgeschaut, es 
gibt doch eine Möglichkeit über das Polynom: Das Ganze nennt sich 
Parameterschätzverfahren von ARX-Modellen.

Hierzu muss Du eine Messmatrix erstellen und die Pseudoinverse daraus, 
dann bekommst Du anscheinend die Parameter des Polynoms raus, die sollte 
auch mit Störgrößen funktionieren.

Der Aufwand liegt eigentlich darin die Pseudo-Inverse zu berechnen.

Gruß
Volly

[Chris (Gast)]

Danke für die Antworten

die Grundstrucktur der Regelstrecke ist bekannte. Es geht wirklich darum 
die Parameter zu bestimmen.

Das Beispiel aus dem Handbuch habe ich auch schon in Matlab ausprobiert. 
Es nutzt ja auch die Methode der kleinsten Quadrate. Es hat sich 
gezeigt, dass es nicht funktieniert wenn es z.b. mit einem Sinus statt 
Random Werten angeregt wird.

Das Java Skript sieht schon gut aus ich versuche das beschriebene um zu 
setzen.

[volly (Gast)]

Hallo Chris,

ich habe mal die ARX auf Scilab getestet, und eigentlich funktioniert es 
auch mit dem Sinus - es darf allerdings kein reiner Sinus sein, sondern 
muss ein bisschen Rauschen haben.

Das Problem liegt wahrscheinlich in dem Lösen des Gleichungssystem mit 
der Pseudoinversen Matrix, diesbezüglich existieren noch andere 
Verfahren, die ich noch nicht geprüft habe.

Anbei mal der Scilab-Code (fast so wie Matlab).

Gruß
Volly

--

function u=ut(t), u=sin(50*t)+0.00001*rand();endfunction

s=poly(0,'s');                      // definiert s als Polynomvariable

Z = 30-0.47407*s-27.36148*s*s;
N = 1-0.3155*s+5.44933*s*s;

t=0:0.001:0.2;
e=[];
for i=1:size(t,2)
  e(i)=ut(t(i));
end

Gs=syslin('c',Z,N)
yk=csim(e',t,Gs);

N = size(t,2)-2;
n = 2;
m = 2;
nmp1 = n+m+1;
M = zeros(N,nmp1);
y = zeros(N,1);

for k1 = 1 : N
  y(k1) = yk(k1+2);
end

for k1 = 1 : N
  M(k1,1) = -yk(k1+1);
  M(k1,2) = -yk(k1);
  M(k1,3) = e(k1+2);
  M(k1,4) = e(k1+1);
  M(k1,5) = e(k1);
end

p = pinv(M) * y

Z = p(3)+p(4)*s+p(5)*s*s;
N = 1+p(1)*s+p(2)*s*s;

Gs=syslin('c',Z,N)
yk1=csim(e',t,Gs);

xbasc(0);
subplot(3,1,1)
xtitle("Eingangsgröße")
plot2d(t,[e],2)
subplot(3,1,2)
xtitle("Gegebene Ausgangsgröße")
plot2d(t,[yk'],3)
subplot(3,1,3)
xtitle("LS Ausgangsgröße")
plot2d(t,[yk1'],5)

[volly (Gast)]

Hallo Chris,

Du hattest doch Recht. Da die Parameter doch irgendwie unterschiedlich 
waren, habe ich mir mal die Bode-Diagramme berechnet - und die sind 
vollständig unterschiedlich.

Gruß
Volly

[Chris (Gast)]

Hallo Volly,

ich werde gleich mal eine Matlab funktion ausprobieren die ich hier 
gefunden habe:

Beitrag "Matlab Sprungantwort Fitting"

ich habe heute folgendes in Octave ausprobiert (werde es auch noch in 
Matlab testen):
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
%Ident TEST Octave
clear all

T = 0.1;
t = [ 0: T :10];
f = 1;
uk = sin(2*pi*f*t);
%uk = rand(1,length(t))*10;      % Eingangssignal
t = t';
uk = uk';

rk = rand(1,length(t))*0.5;      % Störsignal (auskommentier)

K = 5;
T1 = 2;
num = [K];    % b
den = [T1,1];  %a0 a1
GS = tf(num, den);
%sysout(GS);

yk = lsim(GS,uk,t) ;%+ rk      % Ausgangssignal
uk = uk';
figure(1)
subplot(211), plot(t,uk);
subplot(212), plot(t,yk);

N = length(uk) - 1
m = 1;            % Ordnung des Systems
d = 0;            % Totzeit
psi = zeros (N, m);       % Meßmatrix
y = zeros (N, 1);         % Ausgangsgröße
for i = 1:N
    y(i) = yk(m+d+i);    %Belegung der Ausgangsgröße
end

for i = 1:N             %Belegung der Meßmatrix
    %psi(i,1) = -yk(m+d+i-1);
    psi(i,1) = -yk(d+i);      % a Nenner
    %psi(i,3) =  uk(m+i-1);
    psi(i,2) =  uk(i);        % b Zähler
end

omega = pinv(psi) * y
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Das hat aber auch nur mit ein paar Werten funktioniert.
Ich habe mir das aus dem Taschenbuch Beispiel, und dem Buch 
"Identifikation Dynamischer Systeme" von Isermann zusammen gereimt.

Melde mich morgen mit neuen Ergebnissen :-)

Gruß

Chris