|
Autor |
Adams-Moulton-Verfahren |
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Themenstart: 2021-11-30
|
https://www.matheplanet.com/matheplanet/nuke/html/uploads/b/54010_Unbenltkgjt_lgtannt.PNG
Hallo zusammen, könnte mir jemand bitte dabei helfen?
Ich muss diese Formel auf Matlab implementieren :
https://www.matheplanet.com/matheplanet/nuke/html/uploads/b/54010_Unbentogijtoigjtoannt.PNG
Und das ist die Tabelle für beta und k = 1,2,3 :
https://www.matheplanet.com/matheplanet/nuke/html/uploads/b/54010_Unbekhglgtjltgtnannt.PNG
Mein Problem ist dass das Verfahren implizit ist und weiß nicht was ich statt f(x_1,y_1) einsetzen soll
Danke im Voraus!
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.1, vom Themenstarter, eingetragen 2021-11-30
|
https://www.matheplanet.com/matheplanet/nuke/html/uploads/b/54010_71472e10-1008-477c-b9e6-abd4f7da3e1f.jpg
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.2, eingetragen 2021-11-30
|
\quoteon(2021-11-30 20:32 - s-amalgh im Themenstart)
Mein Problem ist dass das Verfahren implizit ist und weiß nicht was ich statt f(x_1,y_1) einsetzen soll
\quoteoff
Du musst ausnutzen, dass hier $f(x_i,y_i)=A\,y_i$ linear in $y_i$ ist.
Wenn du das in die Gleichung, die $y_{i+1}$ implizit definiert, einsetzt, erhältst du ein lineares Gleichungssystem für $y_{i+1}$.
--zippy
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.3, vom Themenstarter, eingetragen 2021-11-30
|
Danke erstmal für deine Antwort!
Wie kann ich das in Matlab implementieren ?
Könntest du mir bitte bei dem ersten Schritt helfen?
Dnke im Voraus!
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.4, eingetragen 2021-11-30
|
\quoteon(2021-11-30 22:37 - s-amalgh in Beitrag No. 3)
Wie kann ich das in Matlab implementieren ?
\quoteoff
Matlab kann lineare Gleichungssysteme lösen, das hast du schonmal verwendet:
\quoteon(2021-11-06 16:25 - s-amalgh in Beitrag No. 24)
\sourceon Matlab
z = -m\z
\sourceoff
\quoteoff
Hier lautet das Gleichungssystem $\left(1-\frac h2A\right)y_1 =
\left(1+\frac h2A\right)y_0$, wobei $1$ die $2\times2$-Einheitsmatrix bezeichnet.
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.5, vom Themenstarter, eingetragen 2021-11-30
|
wie kommst du auf (1-h/2*A) ?
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.6, eingetragen 2021-11-30
|
\quoteon(2021-11-30 22:50 - s-amalgh in Beitrag No. 5)
wie kommst du auf (1-h/2*A) ?
\quoteoff
Die ursrpüngliche Gleichung lautet wegen $\beta_0^*=\beta_1^*=\frac12$$$
y_1 = y_0 + \frac h2\Bigl[f(x_1,y_1)+f(x_0,y_0)\Bigr] =
y_0 + \frac h2\Bigl[A\,y_1+A\,y_0\Bigr] \;,
$$und Sortieren nach $y_1$ und $y_0$ liefert$$
y_1-\frac h2A\,y_1 = y_0+\frac h2A\,y_0 \quad\iff\quad
\left(1-\frac h2A\right)y_1 = \left(1+\frac h2A\right)y_0 \;. $$
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.7, vom Themenstarter, eingetragen 2021-11-30
|
Achso ok..
könntest du mir bitte zeigen wie ich das in matlab implementieren kann?
nur den ersten Schritt, ich ach dann der Rest allein
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.8, eingetragen 2021-11-30
|
\quoteon(2021-11-30 22:37 - s-amalgh in Beitrag No. 3)
Wie kann ich das in Matlab implementieren ?
\quoteoff
Unter der Annahme, dass a die Matrix $A$, y den Vektor $y_0$ und h die Schrittweite $h$ enthält, ergibt sich $y_1$ als:
\sourceon matlab
(eye(2)-.5*h*a)\((eye(2)+.5*h*a)*y)
\sourceoff
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.9, vom Themenstarter, eingetragen 2021-12-01
|
\sourceon Matlab
function t1
y0 = [1; 0];
a = 0;
b = 10;
h = [1e-1,1e-2,1e-3];
A = [0, 1; -4*pi^2, 0];
f = @(x,y) (A*y);
for k = 1:3
for j = 1:length(h)
[~,y] = adammoult(a, b,y0,h(j),1,f);
end
end
function [x,y] = adammoult(a,b,y0,h,k,rhs)
beta = [1/2 1/2 0 0; ... % k = 1
5/12 8/12 -1/12 0; ... % k = 2
9/24 19/24 -5/24 1/24]; ... % k = 3
x = a:h:b;
N = length(x);
n = length(y0);
y = [y0(:), zeros(n,N-1)];
if k == 1
for i = 1:N-1
y(:,i+1) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,i));
end
elseif k ==2
for i = 2:N-1
y(:,2) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,1));
% Hier überlege ich noch was ich schreiben soll %
end
end
end
\sourceoff
So ? noch nicht fertig aber als idee
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.10, eingetragen 2021-12-01
|
Eine Fallunterscheidung nach dem Wert von $k$ brauchst du eigentlich nicht. Du kannst für alle $k$:
1. Die linke Seite von ...\... aus $\beta_0^*$ sofort berechnen.
2. Die rechte Seite von ...\... in einer Schleife als Summe berechnen.
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.11, vom Themenstarter, eingetragen 2021-12-01
|
wie meinst du das ?
Mein bisheriger Code ist so :
\sourceon Matlab
function t1
y0 = [1; 0];
a = 0;
b = 10;
h = [1e-1,1e-2,1e-3];
A = [0, 1; -4*pi^2, 0];
for k = 1:3
for j = 1:length(h)
[~,y] = adammoult(a, b,y0,h(j),k,A);
end
end
function [x,y] = adammoult(a,b,y0,h,k,A)
beta = [1/2 1/2 0 0; ... % k = 1
5/12 8/12 -1/12 0; ... % k = 2
9/24 19/24 -5/24 1/24]; ... % k = 3
x = a:h:b;
N = length(x);
n = length(y0);
y = [y0(:), zeros(n,N-1)];
if k == 1
for i = 1:N-1
y(:,i+1) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,i));
end
elseif k == 2
y(:,2) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,1));
for i = 2:N-1
y(:,i+1) = (eye(2)-5/12*h*A)\((eye(2)+8/12*h*A)*y(:,i) - 1/12*h*A*y(:,i-1));
end
elseif k == 3
y(:,2) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,1));
y(:,3) = (eye(2)-5/12*h*A)\((eye(2)+8/12*h*A)*y(:,2) - 1/12*h*A*y(:,1));
for i = 3:N-1
y(:,i+1) = (eye(2)-9/24*h*A)\((eye(2)+19/24*h*A)*y(:,i) - 5/24*h*A*y(:,i-1) + 1/24*h*A*y(:,i-2));
end
end
end
end
\sourceoff
|
Profil
|
zippy
Senior  Dabei seit: 24.10.2018 Mitteilungen: 4607
 | Beitrag No.12, eingetragen 2021-12-01
|
Das sieht aus, als könnte es funktionieren (ich habe es mit aber nicht im Detail angeschaut).
Allerdings programmierst du einiges aus, was du durch Schleifen oder Funktionen ersetzen könntest. Ein Symptom dafür ist, dass du die Matrix beta gar nicht verwendest.
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.13, vom Themenstarter, eingetragen 2021-12-01
|
Könntst du bitte mein Code bearbeiten?
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.14, vom Themenstarter, eingetragen 2021-12-01
|
wie kann ich jetzt plotten?
plot(h,y(1,:)) so?
wenn ich run mache zeigt dass die beide Vektoren die gleiche Länge sein müssen. Wie mach ich das?
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.15, vom Themenstarter, eingetragen 2021-12-01
|
Profil
|
rlk
Senior  Dabei seit: 16.03.2007 Mitteilungen: 11600
Wohnort: Wien
 | Beitrag No.16, eingetragen 2021-12-01
|
Hallo s-amalgh,
überlege Dir, welche Bedeutung $h$ und der derzeit nicht verwendete Vektor x haben.
Servus,
Roland
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.17, vom Themenstarter, eingetragen 2021-12-01
|
Ich habe nicht verstanden was du gemeint hast..
Kannst du mir bitte das mehr erklären?
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.18, vom Themenstarter, eingetragen 2021-12-01
|
Meinst du das so?
für k = 1:
\sourceon Matlab
h = [1e-1,1e-2,1e-3];
#for k = 1:3
for j = 1:length(h)
[~,y] = adammoult(a, b,y0,h(j),1,A)
plot(0:h(j):10,y(1,:))
hold on
end
\sourceoff
|
Profil
|
rlk
Senior  Dabei seit: 16.03.2007 Mitteilungen: 11600
Wohnort: Wien
 | Beitrag No.19, eingetragen 2021-12-01
|
Hallo s-amalgh,
die Funktion plot() erwartet zwei Vektoren mit den x- und y-Koordinaten der zu zeichnenden Punkte. Die Zahl $h$ ist die Schrittweite, also der horizontale Abstand zweier benachbarter Punkte.
Die Funktion adammoult() liefert zwei Ergebnisse zurück, von denen Du das erste nicht verwendest.
Ich hoffe, das hilft Dir,
Roland
[Die Antwort wurde nach Beitrag No.17 begonnen.]
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.20, vom Themenstarter, eingetragen 2021-12-01
|
für k = 1:
\sourceon Matlab
h = [1e-1,1e-2,1e-3];
#for k = 1:3
for j = 1:length(h)
[~,y] = adammoult(a, b,y0,h(j),1,A)
plot(0:h(j):10,y(1,:))
hold on
end
\sourceoff
Meinst du das so oder wie?
|
Profil
|
rlk
Senior  Dabei seit: 16.03.2007 Mitteilungen: 11600
Wohnort: Wien
 | Beitrag No.21, eingetragen 2021-12-01
|
Hallo s-amalgh,
nein, ich meinte, dass Du den Vektor x , den Du hoffentlich nicht ohne Grund in der Funktion adammoult berechnen lässt, verwendest wofür er gedacht ist.
Servus,
Roland
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.22, vom Themenstarter, eingetragen 2021-12-01
|
Du meintest das so oder:
Für k = 1:
\sourceon Matlab
h = [1e-1,1e-2,1e-3];
for j = 1:length(h)
[x,y] = adammoult(a, b,y0,h(j),1,A)
plot(x,y(1,:))
hold on
end
\sourceoff
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.23, vom Themenstarter, eingetragen 2021-12-01
|
Mein Code :
\sourceon Matlab
function Alghabra07
function [x,y] = adammoult(a,b,y0,h,k,A)
% Gitterpunkte der N"aherungen
x = a:h:b;
% Anzahl der Gitterpunkte
N = length(x);
% Dimension der Dgl
n = length(y0);
% Initialisierung der N"aherungen
y = [y0(:), zeros(n,N-1)];
%Adams-Moulton-Verfahren
% K = 1
if k == 1
for i = 1:N-1
y(:,i+1) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,i));
end
% K = 2
elseif k == 2
y(:,2) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,1));
for i = 2:N-1
y(:,i+1) = (eye(2)-5/12*h*A)\((eye(2)+8/12*h*A)*y(:,i) - 1/12*h*A*y(:,i-1));
end
% K = 3
elseif k == 3
y(:,2) = (eye(2)-.5*h*A)\((eye(2)+.5*h*A)*y(:,1));
y(:,3) = (eye(2)-5/12*h*A)\((eye(2)+8/12*h*A)*y(:,2) - 1/12*h*A*y(:,1));
for i = 3:N-1
y(:,i+1) = (eye(2)-9/24*h*A)\((eye(2)+19/24*h*A)*y(:,i) - 5/24*h*A*y(:,i-1) + 1/24*h*A*y(:,i-2));
end
end
end
% Beispiel
y0 = [1; 0];
% x_0
a = 0;
% x_end
b = 10;
% Schrittweiten
h = [1e-1,1e-2,1e-3];
% Matrix
A = [0, 1; -4*pi^2, 0];
for k = 1:3
for j = 1:length(h)
[x,y] = adammoult(a, b,y0,h(j),k,A)
figure(k);
plot(x,y(1,:));
set(gca,'FontSize',24);
title(['k = ',num2str(k)]);
hold on;
end
hold off;
legend({'h(1)','h(2)','h(3)'},'Location','southwest');
legend('boxoff');
end
end
\sourceoff
Und das sind die Plots :
https://www.matheplanet.com/matheplanet/nuke/html/uploads/b/54010_Unbenkflghlbvghlitbtgrtannt.PNG
Habe ich irgendein Fehler?
|
Profil
|
rlk
Senior  Dabei seit: 16.03.2007 Mitteilungen: 11600
Wohnort: Wien
 | Beitrag No.24, eingetragen 2021-12-05
|
Hallo s-amalgh,
ich habe mir Dein Programm nicht im Detail angesehen, aber die Plots sehen so aus als wenn es funktionieren würde. In der Legende würde ich die Werte der verwendeten Schrittweite statt der wenig informativen Zeichenketten "h(1)" usw. ausgeben.
Servus,
Roland
|
Profil
|
s-amalgh
Wenig Aktiv  Dabei seit: 16.12.2020 Mitteilungen: 374
 | Beitrag No.25, vom Themenstarter, eingetragen 2021-12-05
|
Danke :)
Könntest du mir bitte dabei helfen?
https://www.matheplanet.com/matheplanet/nuke/html/viewtopic.php?topic=256608&start=0&lps=1863759#v1863759
Es wäre sehr nett von dir, da ich das heute abgeben muss..
Danke im Voraus!
|
Profil
|
s-amalgh hat die Antworten auf ihre/seine Frage gesehen. | s-amalgh wird per Mail über neue Antworten informiert. |
|
All logos and trademarks in this site are property of their respective owner. The comments are property of their posters, all the rest © 2001-2023 by Matroids Matheplanet
This web site was originally made with PHP-Nuke, a former web portal system written in PHP that seems no longer to be maintained nor supported. PHP-Nuke is Free Software released under the GNU/GPL license.
Ich distanziere mich von rechtswidrigen oder anstößigen Inhalten, die sich trotz aufmerksamer Prüfung hinter hier verwendeten Links verbergen mögen. Lesen Sie die
Nutzungsbedingungen,
die Distanzierung,
die Datenschutzerklärung und das Impressum.
[Seitenanfang]
|