Iteration < Interpol.+Approx. < Numerik < Hochschule < Mathe < Vorhilfe
|
| Status: |
(Frage) beantwortet  | | Datum: | 18:02 Mi 11.02.2009 | | Autor: | Tbasket |
Hallo zusammen,
ich muss folgende implizit gegebene Funktion per iteration lösen.
P(X<=d-Q,Y>Q)=c
X und Y sind normalverteilte Zufallsvariablen, d ist konstant und gegeben und c auch. Nun ist die Frage wie ich Q per Iteration herausbekomme. Kann man hier Newton verwenden? Die Wahrscheinlichkeiten könnte ich ja als Doppelintegral schreiben..
Ich habe diese Frage in keinem Forum auf anderen Internetseiten gestellt
Wäre sehr dankbar für Eure Hilfe!
Dateianhänge:
Anhang Nr. 1 (Typ: pdf) [nicht öffentlich]
Anhang Nr. 2 (Typ: pdf) [nicht öffentlich]
|
|
| |
|
> Hallo zusammen,
>
> ich muss folgende implizit gegebene Funktion per Iteration
> lösen.
> P(X<=d-Q,Y>Q)=c
>
> X und Y sind normalverteilte Zufallsvariablen, d ist
> konstant und gegeben und c auch. Nun ist die Frage wie ich
> Q per Iteration herausbekomme. Kann man hier Newton
> verwenden? Die Wahrscheinlichkeiten könnte ich ja als
> Doppelintegral schreiben..
Wahrscheinlich ist noch angegeben, dass X und Y
unabhängig sind. Unter dieser Annahme hätte man
also die Gleichung:
[mm] $\underbrace{\integral_{-\infty}^{d-Q}f_X\ dx\ *\integral_{Q}^{\infty}f_Y\ dy\ -\ c}_{N(Q)}\ [/mm] =\ 0$
Um nach Newton eine Nullstelle der Funktion N zu
bestimmen, braucht man ihre Ableitung:
$\ [mm] N'(Q)=\bruch{d}{dQ}\left(\integral_{-\infty}^{d-Q}f_X\ dx\right)\ *\integral_{Q}^{\infty}f_Y\ [/mm] dy\ +\ [mm] \integral_{-\infty}^{d-Q}f_X\ [/mm] dx\ [mm] *\bruch{d}{dQ}\left(\ \integral_{Q}^{\infty}f_Y\ dy\right)$
[/mm]
$\ =\ [mm] -f_X(d-Q)*\integral_{Q}^{\infty}f_Y\ [/mm] dy\ +\ [mm] \integral_{-\infty}^{d-Q}f_X\ [/mm] dx\ [mm] *(-f_Y(Q))$
[/mm]
$\ =\ [mm] -f_X(d-Q)*\integral_{Q}^{\infty}f_Y\ [/mm] dy\ -\ [mm] f_Y(Q)*\integral_{-\infty}^{d-Q}f_X\ [/mm] dx\ $
Die Integrationen für $\ N(Q)$ und $\ N'(Q)$ könnte man
natürlich noch durch die entsprechenden Ausdrücke
mit den Verteilungsfunktionen ersetzen, falls für
diese vielleicht eine Approximationsformel genutzt
werden soll.
Und dann geht der übliche Newton-Algorithmus los:
Startwert [mm] Q_0 [/mm] wählen und mittels [mm] Q_{k+1}=Q_k-\bruch{N(Q_k)}{N'(Q_k)}
[/mm]
weitere Werte [mm] Q_1, Q_2, Q_3 [/mm] ... berechnen und ihre
Konvergenz beobachten.
Das Ganze könnte recht aufwendig werden, aber ich
darf wohl annehmen, dass das nicht für Handrechnung
gedacht ist ... 
Gruß Al-Chwarizmi
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 09:45 Do 12.02.2009 | | Autor: | Tbasket |
Vielen vielen Dank erstmal für die Antwort, wirklich super das Ihr weiterhelft.
Leider habe ich vergessen zu erwähnen, dass die bivariaten Zufallsvariablen unabhängig. Weißt Du trotzdem wie ich vorgehen könnte?
Das Problem will ich als Iteration in Maple lösen, also nicht per Handlosung. Wenn ich die Wahrscheinlichkeiten ausschreibe habe ich allerdings ein Doppelintegral und ich weiß nicht ob das über Newton noch möglich ist. Oder doch, weil ich mit Q nur eine zu suchende Variable habe?
Vielen Dank im Voraus
|
|
|
| |
|
> Leider habe ich vergessen zu erwähnen, dass die bivariaten
> Zufallsvariablen unabhängig. Weißt Du trotzdem wie ich
> vorgehen könnte?
was heißt da "trotzdem": gerade weil sie unabhängig sind,
ist das Ganze viel einfacher als andernfalls bzw. überhaupt
machbar !
> Das Problem will ich als Iteration in Maple lösen, also
> nicht per Handlösung. Wenn ich die Wahrscheinlichkeiten
> ausschreibe habe ich allerdings ein Doppelintegral und ich
> weiß nicht ob das über Newton noch möglich ist. Oder doch,
> weil ich mit Q nur eine zu suchende Variable habe?
Hallo Mario,
ich denke, ich habe schon ziemlich klar notiert, wie
du vorgehen könntest. Der wichtigste Punkt war, dass
man die Ableitung des in zwei Faktorintegrale zerlegbaren
Doppelintegrals via Produktregel leicht hinschreiben kann.
Die Ableitung N', die man für Newton braucht, steht also
zur Verfügung. Nur kommen darin wie in N selber noch
Normalverteilungs-Integrale vor. Ich weiss nicht, ob die
in Maple schon pfannenfertig vorliegen - aber das sollte
wohl zu machen sein. Was du brauchst, ist im Prinzip
nur die Standard-Normalverteilungs-Funktion
[mm] $\Phi(z):=\bruch{1}{\wurzel{2\,\pi}}\integral_{-\infty}^{z}e^{-\bruch{x^2}{2}}\ [/mm] dx$
von der aus man die [mm] \mu-\sigma- [/mm] Normalverteilung $\ F(x)$ erhält
durch:
[mm] F(\mu,\sigma,x):=\Phi\left(\bruch{x-\mu}{\sigma}\right)
[/mm]
Es lohnt sich wohl auch, die Dichtefunktion
$\ [mm] f(\mu,\sigma,x):=\bruch{1}{\sigma*\wurzel{2\,\pi}}*exp\left(-\bruch{1}{2}\left(\bruch{x-\mu}{\sigma}\right) ^2\right)$
[/mm]
bereitzustellen. Mit diesen Funktionen lassen sich die
Formeln für N(Q) und N'(Q) darstellen.
Übrigens: ich bin einfach deiner Idee nachgegangen,
den Tipp "Iteration" als Newton-Iteration zu inter-
pretieren. Ich weiß nicht, ob der Aufgabensteller das
auch so gemeint hat - aber wir sehen jetzt, dass es
offenbar möglich ist, obwohl nicht so ganz einfach.
LG Al-Chw.
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 12:11 Do 12.02.2009 | | Autor: | Tbasket |
Oh nein, ich war wohl nicht richtig ausgeschlafen als ich die ANtwort geschrieben habe. Ich hoffe, DU hast nochmals die Geduld zu antworten. Die Variablen sind NICHT unabhängig, dadurch verkompliziert sich die AUfageb sicher erheblich.....
Newton ist nicht vorgegeben, es steht mir frei wie ich die Iteration mache, wenn es eine einfachere gibt wäre ich dir sehr dankbar wenn du mir sagen könntest wie. EIne Ableitung von einem Doppelintegral in MAple habe ich nicht hinbekommen....
Vielen vielen Dank im Voraus> > Leider habe ich vergessen zu erwähnen, dass die bivariaten
> > Zufallsvariablen unabhängig. Weißt Du trotzdem wie ich
> > vorgehen könnte?
>
> was heißt da "trotzdem": gerade weil sie unabhängig sind,
> ist das Ganze viel einfacher als andernfalls bzw.
> überhaupt
> machbar !
>
> > Das Problem will ich als Iteration in Maple lösen, also
> > nicht per Handlösung. Wenn ich die Wahrscheinlichkeiten
> > ausschreibe habe ich allerdings ein Doppelintegral und ich
> > weiß nicht ob das über Newton noch möglich ist. Oder doch,
> > weil ich mit Q nur eine zu suchende Variable habe?
>
>
> Hallo Mario,
>
> ich denke, ich habe schon ziemlich klar notiert, wie
> du vorgehen könntest. Der wichtigste Punkt war, dass
> man die Ableitung des in zwei Faktorintegrale zerlegbaren
> Doppelintegrals via Produktregel leicht hinschreiben
> kann.
> Die Ableitung N', die man für Newton braucht, steht also
> zur Verfügung. Nur kommen darin wie in N selber noch
> Normalverteilungs-Integrale vor. Ich weiss nicht, ob die
> in Maple schon pfannenfertig vorliegen - aber das sollte
> wohl zu machen sein. Was du brauchst, ist im Prinzip
> nur die Standard-Normalverteilungs-Funktion
>
> [mm]\Phi(z):=\bruch{1}{\wurzel{2\,\pi}}\integral_{-\infty}^{z}e^{-\bruch{x^2}{2}}\ dx[/mm]
>
> von der aus man die [mm]\mu-\sigma-[/mm] Normalverteilung [mm]\ F(x)[/mm]
> erhält
> durch:
>
> [mm]F(x):=\Phi\left(\bruch{x-\mu}{\sigma}\right)[/mm]
>
> Es lohnt sich wohl auch, die Dichtefunktion
>
> [mm]\ f(\mu,\sigma,x):=\bruch{1}{\sigma*\wurzel{2\,\pi}}*exp\left(-\bruch{1}{2}\left(\bruch{x-\mu}{\sigma}\right) ^2\right)[/mm]
>
> bereitzustellen. Mit diesen Funktionen lassen sich die
> Formeln für N(Q) und N'(Q) darstellen.
>
> Übrigens: ich bin einfach deiner Idee nachgegangen,
> den Tipp "Iteration" als Newton-Iteration zu inter-
> pretieren. Ich weiß nicht, ob der Aufgabensteller das
> auch so gemeint hat - aber wir sehen jetzt, dass es
> offenbar möglich ist, obwohl nicht so ganz einfach.
>
> LG Al-Chw.
|
|
|
| |
|
> Die Variablen sind NICHT unabhängig, dadurch verkompliziert
> sich die Aufgabe sicher erheblich.....
Allerdings. Aber was weisst du denn noch über X und Y ?
Und: kannst du etwas über den Hintergrund der Aufgabe
verraten ?
Mit den unabhängigen normalverteilten X und Y wäre es
ja jetzt eigentlich ganz nett geworden: die Aufgabe hat
mich inzwischen gepackt und ich möchte sie zu Ende bringen.
> Newton ist nicht vorgegeben, es steht mir frei wie ich die
> Iteration mache, wenn es eine einfachere gibt wäre ich dir
> sehr dankbar wenn du mir sagen könntest wie.
Da sehe ich im Moment nichts Passendes, ausser viel-
leicht schrittweise Änderungen von Q in der Art von
 Regula falsi.
LG Al-Chw.
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 14:22 Do 12.02.2009 | | Autor: | Tbasket |
Toll, das Du hilfst!!! Und sorry für meine verspätete Antwort
Also die Zufallsvariablen sind normalverteilt (mu und varianz und Korrelation sind gegeben). Die Aufgabe ist aus dem Bereich Produktion. d ist die Nachfrage und Q die optimale Produktionsmenge, c die Kosten einer zusätzlichen Einheit. Die optimale Menge ist nur explizit gegeben und hängt von den bivariaten Zufallsvariablen X und Y ab.
Wie ich Iteriere ist mir völlig überlassen. Es sollte nur in Maple machbar sein. Newton geht hier wahrscheinlich nicht? Deine vorgeschlagene Methode schaue ich mir direkt mal an, ich denke sie ist leichter, da hier keine ABleitungen vorkommen. (Ich bräuchte ein Verfahren für das es zum einen kein Problem ist, dass eine bivariate W. vorkommt bzw. das Doppelintegral und zum anderen das die Iteration bei null anfangen muss, da man keine Ahnung hat wo Q liegt um es näher einzuschränken....)
Wenn Du weitere Details brauchst einfach fragen
> > Die Variablen sind NICHT unabhängig, dadurch verkompliziert
> > sich die Aufgabe sicher erheblich.....
>
> Allerdings. Aber was weisst du denn noch über X und Y ?
>
> Und: kannst du etwas über den Hintergrund der Aufgabe
> verraten ?
>
> Mit den unabhängigen normalverteilten X und Y wäre es
> ja jetzt eigentlich ganz nett geworden: die Aufgabe hat
> mich inzwischen gepackt und ich möchte sie zu Ende
> bringen.
>
> > Newton ist nicht vorgegeben, es steht mir frei wie ich die
> > Iteration mache, wenn es eine einfachere gibt wäre ich dir
> > sehr dankbar wenn du mir sagen könntest wie.
>
> Da sehe ich im Moment nichts Passendes, ausser viel-
> leicht schrittweise Änderungen von Q in der Art von
> Regula falsi.
>
>
> LG Al-Chw.
|
|
|
| |
|
Hallo Mario,
danke für die Zusatzinformationen. Eben habe ich
übrigens die "alte" Version mit den unabhängigen
X und Y in Mathematica zum Laufen gebracht.
Funktioniert offenbar prima. Um die Formeln etwas
handlicher zu machen, habe ich noch ein paar
Abkürzungen J,K,L,M eingeführt, und statt N(Q)
und N'(Q) heisst es jetzt P[Q] und P1[Q].
Hier das Notebook als Text:
![[a]](/images/paperclip.gif "Dateianhänge") Datei-Anhang
Die Sache mit den korrelierten Variablen schaue ich
mir einmal an, garantiere aber nichts ...
Gruß Al-Chw.
Dateianhänge:
Anhang Nr. 1 (Typ: m) [nicht öffentlich]
Anhang Nr. 2 (Typ: nb) [nicht öffentlich]
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 16:05 Do 12.02.2009 | | Autor: | Tbasket |
Vielen Dank! Echt toll das du es ausprobiert hast
Ich habe leider kein Mathematika und kenne Tex leider nicht, dazu kommt das ich noch einen Apple rechner habe.... Konnte aber den Code einsehen als html seite einsehen
Hast du Newton benutzt oder ein anderes Iterationsverfahren ?
Wenn du es für den korrelierten Fall hinbekommen würdest wäre mein Problem gelöst und würde nicht mehr an dieser STelle hängenbleiben, wie die letzten 2 Wochen :(
|
|
|
| |
|
> Ich habe leider kein Mathematica und kenne TeX leider
> nicht,
Ich habe das bereits auf gewöhnliches Textformat
gebracht und neu reingestellt.
> dazu kommt das ich noch einen Apple rechner habe....
ich auch - und ich werde dabei bleiben
(und insbesondere entschuldige ich mich nicht
dafür, ich bin ja schließlich nicht blöd !)
> Hast du Newton benutzt oder ein anderes
> Iterationsverfahren ?
Newton, genau so wie vorbereitet.
> Wenn du es für den korrelierten Fall hinbekommen würdest
> wäre mein Problem gelöst
Wenn nur das Doppelintegral berechnet werden
kann, sollte das mit Regula falsi kaum ein Problem
sein. Schick doch bitte einmal noch das Doppelintegral
mit den Bezeichnungen, die du verwendet hast !
(ich habe z.B. die Bezeichnungen mux, sigmax, muy,
sigmay etc. verwendet)
Al
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 16:27 Do 12.02.2009 | | Autor: | Tbasket |
Hier habe ich die Allgemeine Form der bivariate normal densinity function:
#Bivariate normal density
[mm] f:=(x,y)->exp((-1/(2*(1-rho^2)))*(((x-mu[1])/sigma[1])^2-2*rho*(x-mu[1])*(y-mu[2])/(sigma[1]*
[/mm]
[mm] sigma[2])+((y-mu[2])/sigma[2])^2))/(2*Pi*sigma[1]*sigma[2]*sqrt(1-rho^2));
[/mm]
# obere und untere Integrationsgrenze
k:=max(mu[1],mu[2])+3*sigma[1]+3*sigma[2]:
# zweimal die obig definierte funktion integriert
F:=(x,y)->evalf[5](int(int(f(s,t),t=-k..y),s=-k..x)):
Kann man nur im ersten Fall newton benutzen oder auch im zweiten?
Ich versuche gerade noch dein ersten Anhang zu verstehen :)
|
|
|
| |
|
> Hier habe ich die Allgemeine Form der bivariate normal
> density function:
>
> #Bivariate normal density
>
> [mm]f:=(x,y)->exp((-1/(2*(1-rho^2)))*(((x-mu[1])/sigma[1])^2-2*rho*(x-mu[1])*(y-mu[2])/(sigma[1]*[/mm]
>
> [mm]sigma[2])+((y-mu[2])/sigma[2])^2))/(2*Pi*sigma[1]*sigma[2]*sqrt(1-rho^2));[/mm]
O.K., habe das Ding eben auch gerade nachgeschlagen.
>
> # obere und untere Integrationsgrenze
> k:=max(mu[1],mu[2])+3*sigma[1]+3*sigma[2]:
>
> # zweimal die obig definierte funktion integriert
> F:=(x,y)->evalf[5](int(int(f(s,t),t=-k..y),s=-k..x)):
Aha, das ist wohl eine Methode, die Integration
aussen abzuschneiden, damit man keine unendlichen
Integrationsgrenzen braucht ...
> Kann man nur im ersten Fall Newton benutzen oder auch im
> zweiten?
Hier macht die Newton-Iteration keinen Sinn, weil
die Ableitung viel zu kompliziert wird.
Wir haben doch jetzt Folgendes:
Es soll die Gleichung
$\ P(X<d-q\ [mm] \wedge\ [/mm] Y>q)\ =\ c$
nach q aufgelöst werden. Mit dem Doppelintegral
geschrieben heisst das (ich bleibe mal bei den
uneigentlichen Integralen):
$\ [mm] \integral_{-\infty}^{d-q}\left(\ \integral_{q}^{\infty}f(x,y)\,dy\right)dx\ [/mm] -\ c\ =\ 0$
Mit deinem obigen k kannst du die unendlichen
Grenzen ersetzen (falls Maple das nicht von
selbst irgendwie schafft) und hast:
$\ [mm] \underbrace{\integral_{-k}^{d-q}\left(\ \integral_{q}^{k}f(x,y)\,dy\right)dx\ -\ c}_{p(q)}\ \approx\ [/mm] 0$
Wenn diese Funktion p(q) einmal implementiert
ist (am besten testest du sie an ein paar Beispielen,
um sicher zu sein, dass sie tut, was sie soll), fehlt
nur noch der Algorithmus für Regula falsi.
Dazu zuerst noch eine wichtige Frage: Ist p(q)
streng monoton ? Nach meiner Überlegung müsste
es wohl so sein, aber ich bin nicht ganz sicher.
Unter dieser Voraussetzung, und da p stetig ist,
existiert genau eine Nullstelle, und die findet man
mit Regula falsi.
Zum Starten braucht man zwei Anfangswerte
[mm] q_0 [/mm] und [mm] q_1, [/mm] am besten so, dass [mm] p(q_0) [/mm] und [mm] p(q_1)
[/mm]
entgegengesetzte Vorzeichen haben. Vielleicht
empfehlen sich [mm] q_0:=-k [/mm] und [mm] q_1:=k [/mm] !
Hier der Alchwarizmus in Pascal, das lässt sich sicher
leicht in Maple umsetzen:
epsilon:=1E-6;
qold:= [mm] q_0 [/mm] ;
qact:= [mm] q_1 [/mm] ;
while abs(qact-qold)>epsilon do
begin
pold:=p(qold); (alte Werte)
pact:=p(qact); (aktuelle Werte)
qnew:=qact-(qact-qold)/(pact-pold)*p(qact); (neuen Wert berechnen)
qold:=qact;
qact:=qnew;
end;
write(qact,p(qact))
Gruß Al-Chw.
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 18:05 Do 12.02.2009 | | Autor: | Tbasket |
Ganz toll! Ich denke so muss ich es hinbekommen. Ich werde es morgen versuchen für diesen Bespielfall in Maple zu programmieren und dann hier online stellen, es wäre natürlich super wenn du nochmal drüber schauen könntest!
wenn Q steigt steigt auch die Wahrscheinlichkeit, ich denke das heißt auch das es monoton ist oder?
Ich habe zu deiner ersten Lösung eine PDF angehängt mit kurzen Fragen (im letzten EIntrag) wenn du nur kurz drüber schauen könntest wäre das toll!
Zu deinem Code noch eine Frage du schreibst
qnew:=qact-(qact-qold)/(pact-pold)*p(qact); gehört da noch ein p vor das (pact-pold)?
Besten Dank und Herzliche Grüße,
Mario
|
|
|
| |
|
> Ich habe zu deiner ersten Lösung eine PDF angehängt mit
> kurzen Fragen (im letzten EIntrag) wenn du nur kurz drüber
> schauen könntest wäre das toll!
schon erledigt; hat ein bisschen gedauert
> Zu deinem Code noch eine Frage du schreibst
> qnew:=qact-(qact-qold)/(pact-pold)*p(qact); gehört da noch
> ein p vor das (pact-pold)?
nee, genau richtig so !
|
|
|
| |
|
| Status: |
(Mitteilung) Reaktion unnötig  | | Datum: | 17:08 Do 12.02.2009 | | Autor: | Tbasket |
Ich habe leider mit meinem "beschränkten Mathematischen Kenntnissen" noch ein paar Fragen zu deiner ersten Lösung und wäre dir sehr dankbar wenn du kurz drüberschauen könntest.
Ich habe die Fragen in deinen Lösungstext geschrieben und nummeriert dass du einfacher antworten kannst
Dateianhänge:
Anhang Nr. 1 (Typ: pdf) [nicht öffentlich]
|
|
|
| |
|
1.# OK
2.# Phi ist die Standardnormalverteilung mit sigma=1
Integrationsvariable x, Integration von [mm] -\infty [/mm] bis z
z ist einfach eine Variable (andere als x), um Phi darzustellen
3.# F ist die Normalverteilungsfunktion mit mu und sigma, die
man durch diese Transformation aus Phi bekommt
4.# OK
J OK
K OK
Für L und M schau dir die Formel nochmals an,
die nach der Ableitung mit Produktregel entstanden ist.
Die noch vorhandenen Integrale lassen sich durch Phi
bzw. durch F ausdrücken.
... und ich merke gerade, dass ich die Funktion F gar
nicht eingesetzt habe. L und M liessen sich also auch
einfacher schreiben:
L[Q]:=F[mux,sigmax,d-Q]
M[Q]:=1-F[muy,sigmay,Q]
P[Q] ist die Funktion, die ich vorher N(Q) genannt
hatte, also die Funktion, deren Nullstelle wir suchen.
Ich musste sie umtaufen, weil das grosse N in Mathe-
matica reserviert ist.
P1[Q] ist die Ableitungsfunktion von P[Q].
Die Formel newton[Q]:=Q-P[Q]/P1[Q] ist einfach die
Newtonsche Rekursionsformel. Der alte Wert von Q
wird dabei durch den neuen ersetzt.
Nachher folgt das Rechenbeispiel:
Startwert 1 eingegeben und dann mehrmals mit
newton[%] aus dem gerade gehabten Wert den
neuen berechnet. Man sieht, dass das Verfahren
nach wenigen Schritten bei Q=1.82146 stehen
bleibt. Dies wäre also (falls wirklich alles richtig
gelaufen ist) die (bzw. eine) Lösung für die oben
definierten Parameter
mux = 4; sigmax = 2; muy = 3; sigmay = 1; d = 5; c = 0.3
LG
|
|
|
| |
|
| Status: |
(Mitteilung) Reaktion unnötig | | Datum: | 18:38 Do 12.02.2009 | | Autor: | Tbasket |
Abermals danke für deine Unterstützung. Toll wie Du so schnell und gut geantwortet hast. Ich versuche deine Ideen beispeilhaft zu implementieren und stelle sie dann hier rein....
|
|
|
| |
|
> Abermals danke für deine Unterstützung. Toll wie Du so
> schnell und gut geantwortet hast. Ich versuche deine Ideen
> beispielhaft zu implementieren und stelle sie dann hier
> rein....
Das war auch wirklich ein interessanter Dialog, der mir sehr
Spass gemacht hat.
Schönen Abend noch !  Al
|
|
|
| |
|
| Status: |
(Frage) beantwortet | | Datum: | 14:15 Sa 14.02.2009 | | Autor: | Tbasket |
Hallo,
wie versprochen hier meine Implementation in Maple als PDF. Ich habe den Code ein wenig anders geschrieben, probiere ihn aber auch noch aus wie Du es vorgeschlagen hast. Die Ergebnisse sind nur logisch wenn sie für Q eine positive Zahl liefern..... Wenn ich den Parameter c2 variere zwischen 0.38 und 2.58 müsste es eigentlich auch immer positive Q Were ergeben. Leider ist dies nicht der fall. Sobald in der Iteration fa und fb nicht mehr unterscheidlich sind spuckt er negative werte aus für Q.
Hast Du eine Idee woran es hängen könnte?
Dateianhänge:
Anhang Nr. 1 (Typ: pdf) [nicht öffentlich]
Anhang Nr. 2 (Typ: pdf) [nicht öffentlich]
|
|
|
| |
|
> Hallo,
>
> wie versprochen hier meine Implementation in Maple als PDF.
> Ich habe den Code ein wenig anders geschrieben, probiere
> ihn aber auch noch aus wie Du es vorgeschlagen hast. Die
> Ergebnisse sind nur logisch wenn sie für Q eine positive
> Zahl liefern..... Wenn ich den Parameter c2 variere
> zwischen 0.38 und 2.58 müsste es eigentlich auch immer
> positive Q Were ergeben. Leider ist dies nicht der fall.
> Sobald in der Iteration fa und fb nicht mehr
> unterscheidlich sind spuckt er negative werte aus für Q.
> Hast Du eine Idee woran es hängen könnte?
woran genau dies liegt, weiss ich auch nicht
Es gibt einen "kleinen" Fehler, der sich natürlich krass
auswirkt: du hast die Reihenfolge der Integrationen
vertauscht. Die innere (von Q bis +k) geht über y, die
äussere über x. Übrigens konnte ich das ganze jetzt
auch programmieren, nicht in Mathematica, aber in
MuPAD. In der Version, die ich hier anhänge, habe ich
die Daten deines Beispiels benützt. Dazu musste ich
natürlich auch noch deinen Parameter p2 einbringen.
Übrigens war meine früher geäusserte Idee, -k und +k
als Startwerte für die Iteration zu nehmen, wohl nicht
so gut. Das führte (bei mir jedenfalls) zu sehr langen
Rechenzeiten.
LG
Datei-Anhang
Dateianhänge:
Anhang Nr. 1 (Typ: pdf) [nicht öffentlich]
|
|
|
| |
|
| Status: |
(Mitteilung) Reaktion unnötig | | Datum: | 17:53 Sa 14.02.2009 | | Autor: | Tbasket |
Hallo,
wie versprochen hier meine Implementation in Maple als PDF. Ich habe den Code ein wenig anders geschrieben, probiere ihn aber auch noch aus wie Du es vorgeschlagen hast. Die Ergebnisse sind nur logisch wenn sie für Q eine positive Zahl liefern..... Wenn ich den Parameter c2 variere zwischen 0.38 und 2.58 müsste es eigentlich auch immer positive Q Were ergeben. Leider ist dies nicht der fall. Sobald in der Iteration fa und fb nicht mehr unterscheidlich sind spuckt er negative werte aus für Q.
Hast Du eine Idee woran es hängen könnte?
Super vielen Dank. Ich habe für deine Werte ebenefalls 76.88 rausbekommen nachdem ich die Integrationsgrenzen richtig gestellt habe......
Ich mache dann mal weiter, kann leider durchaus sein, dass ich deine Hilfe nochmal in Anspruch nehmen muss....
Ich muss mir auch nochmal den Unterschied zwischen meiner R.F. Iteration und deiner R.F. anschauen, verstehe ich noch nicht ganz....
LG
|
|
|
| |
|
> Wenn ich den Parameter c2 variiere
> zwischen 0.38 und 2.58 müsste es eigentlich auch immer
> positive Q Werte ergeben. Leider ist dies nicht der Fall.
> Sobald in der Iteration fa und fb nicht mehr
> unterscheidlich sind spuckt er negative Werte aus für Q.
Wenn fa=fb (exakt), sollte eigentlich noch Schlimmeres
geschehen: Division durch Null !
Um derartige numerische Probleme zu vermeiden, gäbe
es noch ein "primitiveres", aber stabileres Verfahren
als Regula falsi, nämlich die Nullstellensuche durch
fortgesetzte Intervallhalbierung ...
Da könnte man wohl ruhig auch mit [-k;+k] als
Startintervall beginnen. (***)
Das ginge etwa so:
a:= -k; b:=k;
epsilon:=min(sigmax,sigmay)/100; (***)
while abs(a-b)>epsilon do
m:=(a+b)/2;
if p(a)*p(m)>0 then a:=m else b:=m; end_if (***)
print(float(m))
end_while
(***)Bei dieser Aufgabe kann man davon ausgehen,
dass p(-k)<0 und p(+k)>0 ist. Damit ist die Startbedingung
für das Intervallhalbierungsverfahren erfüllt, und da p stetig
ist, muss (mindestens) eine Nullstelle existieren.
(***)Bemerkung: vorher habe ich mit mikroskopischen
Epsilonwerten gerechnet. Das macht kaum Sinn, wenn die
Zahlenwerte, um die es geht, sehr gross sind. Man kann
also Epsilon proportional zu den Standardabweichungen
oder ev. auch proportional zu k wählen.
(***)Edit: hier hatte ich dummerweise zuerst a und b
verwechselt ...
Gruß Al-Chwarizmi
|
|
|
|
