ich bin noch recht neu in der Arbeit mit Matlab und stehe nach langer Internetrecherche vor einem für mich unlösbaren Problem - das für mich aber eigentlich recht simpel wirkt.
Ich möchte die Summe einer kompontenweisen Multiplikation zweier Vektoren bilden. Die Komponenten des einen Vektors sollen dabei noch Teil einer e-Funktion sein. Die mathematische Formel dazu sieht so aus:
N, n und t sind Konstanten
C ist eine Matrix, von der jeweils jedoch nur die Zeile (N-n) verwendet wird
mueh ist ein Vektor
Ein Aufsummieren mittels for-Schleifen geht nur für kleine Werte für N. Wird N größer, so ergeben sich durch die zahlreichen Rundungen Fehler in einem Ausmaß, das die Anwendung unmöglich macht. Allerdings muss es für meine Anwendung auch möglich sein, größere Werte für N einzusetzen.
Die Umsetzung mit for-Schleifen würde so aussehen:
Code:
for count=1:1:101 for g=0:1:N-n
P(count)=P(count)+(C((N+1)-n,g+1)*exp(-t*mueh(N-g)));
end end
Die Variable Count soll dabei den Vektor P durchlaufend indizieren. Die von der oben stehenden Formel abweichende Indizierung habe ich gewählt, weil bei Matlab die Indizes größer 0 sein müssen.
Hat jemand von euch eine Idee, wie man dies ohne Schleife schreiben kann, sodass möglichst nur einmal ("am Ende") gerundet wird. Ich habe mich auch bereits mit sum, cumsum und sym versucht, bin aber nicht weiter gekommen. Gescheitert bin ich dabei meist an der "doppelten Indizierung" innerhalb der Summe - Fehlermeldung ("Error using * Inner matrix dimensions must agree")
Deine Hoffnung, dass bei einer vektorisierten Methode die Lösung nur einmal gerundet wird, ist leider unbegründet. Auch intern rechnen die Prozessoren nur mit IEEE754 doubles. Ob die Schleife in Matlab oder im C- oder Assembler-Code der Libraries durchgeführt wird, macht keinen Unterschied.
Mit
sym
solltest Du die Summe schon symbolisch lösen können.
vpa
hilft auch, mit Zahlen höherer Genauigkeit zu rechnen.
Was bedeutet "kleinere" und "größere" Werte für N? Manche User finden 100 schon groß, weil sie das nicht mehr per Taschenrechner hinbekommen. Andere arbeiten auf einem Rechen-Cluster und finden alles unter 10^12 Elementen nicht erwähnenswert.
Verfasst am: 12.07.2017, 16:30
Titel: Leider kein Erfolg...
Hallo nochmals,
ich habe mich mal an Jans Lösungen versucht. Jan, du lagst mit deiner Einschätzung völlig richtig, die symbolische Lösung ändert nichts...
Leider hilft auch die Verwendung von XSum und vpa nicht weiter. Das Problem ist wohl, dass die C-Werte bei wachsendem N sehr groß werden. Sie weißen ein alternierendes Vorzeichen auf, was dazu führt, dass bei der Summation dennoch bereits kleine Rundungsungenauigkeiten problematisch sind.
Hier einmal der gesamte Code in symbolischer Schreibweise und für N=20. Die Wrete für mueh müssen so eng beeinander liegen.
Code:
% Vorgabe von N, µ und Betrachtungszeitraum (T)
N=20;
mueh=[0.740.750.760.770.780.790.800.810.820.830.840.850.860.870.880.890.900.910.921];
% Anlegen der Matrix C - dabei beachten: Indexwerte der Matlab-Matrix % sind stets um 1 höher als in der ursprünglichen Formel
C=ones(N,N);
% Berechunung der C-Werte for i=2:1:N
for k=1:1:N
zaehler=1;
nenner=1;
for j=0:1:(i-2)
zaehler=zaehler*mueh(N-j);
end for l=0:1:(i-1) if(l~=(k-1))
nenner=nenner*(mueh(N-l)-mueh(N-k+1));
end end
C(i,k)=zaehler/nenner;
end end
% Anlegen der Matrizen für die symbolische Berechnung
CfF=zeros(N,N);
MfF=zeros(N,N);
P=sym('p',[1 N]);
syms t;
% Berechnung der Wahrscheinlichkeiten als symbolische Formel für jedes n for n=1:1:N
Verfasst am: 12.07.2017, 18:17
Titel: Re: Leider kein Erfolg...
Hallo JogiU,
Mit XSum kann man mit 8 facher Genauigkeit summieren, wobei erst das Endergebnis wieder auf einen Double gerundet wird. Mit VPA kannst Du die verwendete Genauigkeit einstellen.
Wieso helfen die beiden Methoden dann nicht weiter?
Bei der Verwendung von sym und vpa ist es wichtig, dass man nicht erst das Ergebnis als symbolisch oder mit hoher Precision berechnet, sondern dass schon alle Summanden so vorliegen müssen.
Wenn in "f(t)=sum((CfF(n,:) .* exp(-MfF(n,:)*t)));" die rechte Seite schon gerundet wird, bringt es nichts, sie hinterher in die sym,bolische Variable P zu schreiben.
Da ich die Symbolic Toolbox nicht habe, kann ich dazu leider keinen Code vorschlagen.
Es stimmt natürlich: Manche Algorithmen sind numerisch instabil und winzige Rundungsfehler werden dramatisch verstärkt. Das liegt dann am Algorithmus und lässt sich auch prinzipiell nicht mit Rechen-Tricks umgehen. Das ist vergleichbar mit dem mechanischen nd optischen Problem, von der Erde aus einen Laser auf den Mond zu richten und dort stabil nur ein 1cm^2 großen Punkt beleuchten zu wollen. Das geht einfach nicht.
abermals danke für deine Antwort!
Ich habe vpa und XSum folgendermaßen verwendet:
Code:
% Vorgabe von N, µ und Betrachtungszeitraum (T)
N=11;
mueh=[0.830.840.850.860.870.880.890.900.910.921];
T=25;
% Anlegen der Matrix C - dabei beachten: Indexwerte der Matlab-Matrix % sind stets um 1 höher als in der ursprünglichen Formel
C=ones(N,N);
% Berechunung der C-Werte for i=2:1:N
for k=1:1:N
zaehler=1;
nenner=1;
for j=0:1:(i-2)
zaehler=zaehler*mueh(N-j);
end for l=0:1:(i-1) if(l~=(k-1))
nenner=nenner*(mueh(N-l)-mueh(N-k+1));
end end
C(i,k)=zaehler/nenner;
end end
% Anlegen der Wahrscheinlichkeitsmatrix der Zählvariable sowie des % Zeitvektors für den Plot
P=zeros(N+1,101);
count=1;
Zeitvektor=zeros(1,101);
% Berechnung der Wahrscheinlichkeiten in Abhängigkeit von t
for t=0:(T/100):T
for n=1:1:N
for g=0:1:N-n
Vektor(g+1)=vpa((C((N+1)-n,g+1)*exp(-t*mueh(N-g))),40);
end
P(n,count)=XSum(Vektor,[],'Knuth2');
Vektor=Vektor*0;
end
%P(N+1,count)=1-(sum(P((1:1:N),count)));
Zeitvektor(1,count)=t;
count=count+1;
end
Ich befürchte, es ist nicht die optimale Lösung, aber etwas besseres ist mir dabei nicht eingefallen, um die hohe Genauigkeit möglichst von Anfang an drin zu haben. Oder doch total daneben? Ein Veränderung (abgesehen von deutlich längerer Rechenzeit) konnte ich auch für verschiedene Werte von vpa nicht erkennen.
Bei der symbolischen Schreibweise dachte ich, die Rundung der e-Funktion wird erst am Ende vorgenommen. Ich habe es mir aber gerade einmal näher angeschaut und tatsächlich, wird der nicht symbolische Teil der E-Funktion direkt nach vorne gezogen und verrechnet. Wenn ich dich richtig verstehe, müsste ich also Matlab mittels der Symbolic Toolbox dazu bringen, vor der Summation nichts aus der e-Funktion "herauszuziehen"?
Irgendetwas läuft schief. XSum akzeptiert nur numerische Double Arrays als Input. Wenn der Code läuft, ist Vektor also keine symbolische Variable und vpa bringt rein gar nichts.
Nun experimentiere mit sym und vpa, bzw mit XSum (aber nicht mit beidem gleichzeitig, denn das sind unterschiedliche Konzepte, die nicht zusammen passen).
Gruß, Jan
Einstellungen und Berechtigungen
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.
RSS