|
Autor |
cg solver in c++ |
|
elliment
Neu  Dabei seit: 20.12.2019 Mitteilungen: 2
 | Themenstart: 2022-08-07
|
Hallo,
ich habe ein Problem bezüglich eines Conjugate Gradient solvers.
Ich versuche folgendes zu implementieren:
https://matheplanet.de/matheplanet/nuke/html/uploads/b/52376_Bild_2022-08-07_205232011.png
Es soll folgendes lösen:
https://matheplanet.de/matheplanet/nuke/html/uploads/b/52376_Bild_2022-08-07_205811857.png
Mein Code sieht folgendermaßen aus:
\sourceon c++
\numberson
void CG_Solver(Array& Sol, Array& RHS, double& hx, double& hy, int& iters, int nx, int ny) {
// set start timer
auto start = std::chrono::high_resolution_clock::now();
double hxs = hx * hx;
double hys = hy * hy;
Array r;
r.resize(nx, ny);
for (int i = 1; i < Sol.nrows() - 1; i++) {
for (int j = 1; j < Sol.ncolms() - 1; j++) {
r(i, j) = RHS(i,j) - (-1 / hxs * (Sol(i + 1, j) - 2 * Sol(i, j) + Sol(i - 1, j)) - 1 / hys * (Sol(i, j + 1) - 2 * Sol(i, j) + Sol(i, j - 1)));
}
}
double rohold = l2_norm_residual(r, RHS, hx, hy);
double rohnew = 0.0;
Array d;
Array z;
d.resize(nx, ny);
z.resize(nx, ny);
double alpha = 0.0;
double sum = 0.0;
double beta = 0.0;
for (int i = 0; i < Sol.nrows(); i++) {
for (int j = 0; j < Sol.ncolms(); j++) {
d(i, j) = r(i, j);
}
}
for (int k = 0; k < iters; k++) {
for (int i = 1; i < Sol.nrows() - 1; i++) {
for (int j = 1; j < Sol.ncolms() - 1; j++) {
z(i, j) = -1 / hxs * (d(i + 1, j) - 2 * d(i, j) + d(i - 1, j)) - 1 / hys * (d(i, j + 1) - 2 * d(i, j) + d(i, j - 1));
sum += d(i, j) * z(i, j);
}
}
alpha = pow(rohold, 2) / sum;
std::cout << "alpha " << alpha << std::endl;
std::cout << "rohold " << rohold << std::endl;
std::cout << "sum " << sum << std::endl;
for (int i = 1; i < Sol.nrows() - 1; i++) {
for (int j = 1; j < Sol.ncolms() - 1; j++) {
Sol(i, j) = Sol(i, j) + alpha * d(i, j);
r(i, j) = r(i, j) - alpha * z(i, j);
}
}
rohnew = l2_norm_residual(r, RHS, hx, hy);
beta = pow(rohnew, 2) / pow(rohold, 2);
for (int i = 1; i < Sol.nrows() - 1; i++) {
for (int j = 1; j < Sol.ncolms() - 1; j++) {
d(i, j) = r(i, j) + beta * d(i, j);
}
}
}
// set end timer
auto end = std::chrono::high_resolution_clock::now();
// calculate execution time
auto duration = std::chrono::duration_cast(end - start);
std::printf("Time to solution: %.6f seconds\n", (double)duration.count() / 1000000);
}
\sourceoff
Das Problem ist, das meine Werte für u sich nicht verändern bzw. ungefähr null sind und ich kann nicht feststellen wo das Problem liegt. Die Rechte Seite sowie die Grenzen wurden richtig eingestellt. Deswegen, denke ich, dass das Problem am solver liegt, ich kann jedoch den Fehler nicht finden.
Hier ist noch der Gnuplot den ich rauskrieg.
https://matheplanet.de/matheplanet/nuke/html/uploads/b/52376_Bild_2022-08-07_210237585.png
ich hoff mir kann einer helfen und auf den Fehler im Code hinweisen.
Danke im Voraus 🙂
|
Profil
|
ligning
Senior  Dabei seit: 07.12.2014 Mitteilungen: 3555
Wohnort: Berlin
 | Beitrag No.1, eingetragen 2022-08-07
|
Hallo,
man sieht, dass die Variable sum nicht initialisiert ist.
Wenn du möchtest, dass sich das jemand näher ansieht, solltest du den Quelltext als Text posten. Oder soll man das etwa abtippen? Am besten den Quelltext in eine Quelltext-Umgebung posten:
\sourceon C++
std::cout << "Hallo Welt";
\sourceoff
|
Profil
|
Delastelle
Senior  Dabei seit: 17.11.2006 Mitteilungen: 2430
 | Beitrag No.2, eingetragen 2022-08-08
|
Hallo elliment!
Ich weiß nicht ob ich helfen kann!
ich habe vor längerer Zeit mal Innere Dirichlet Probleme im Rechteck gelöst bzw. gelöst bekommen.
Für die mir bekannten Fälle konnte man exakte Lösungen finden.
Dabei kamen sin, cos, sinh, cosh Terme zur Anwendung.
Als ich jetzt mal eine Lösung plotten wollte, bekam ich Probleme mit den z=f(x,y) Werten. Diese sind teilweise viel zu groß. Vielleicht Rundungsfehler bei den sinh und cosh-Termen?
Siehe hier:
https://www.matheplanet.de/matheplanet/nuke/html/viewtopic.php?topic=259707&start=0&lps=1885490#v1885490
Vielleicht kannst Du deinen Löser mal mit verschieden exakten Gittern testen! So n_x = 5, 10, 50, 100 oder ähnlichem... Vielleicht erkennt man dabei bereits etwas!
Es gab auch mal einen Matheplanet-Artikel zum CG-Verfahren -
der könnte eventuell auch helfen:
https://www.matheplanet.de/matheplanet/nuke/html/article.php?sid=1244
Viele Grüße
Ronald
|
Profil
|
majoka
Senior  Dabei seit: 25.02.2014 Mitteilungen: 814
 | Beitrag No.3, eingetragen 2022-08-08
|
Hallo elliment,
Du bist ja mittlerweile der Bitte von ligning nachgekommen:
\quoteon(2022-08-07 21:44 - ligning in Beitrag No. 1)
Wenn du möchtest, dass sich das jemand näher ansieht, solltest du den Quelltext als Text posten.
\quoteoff
Besser wäre es jedoch, wenn Du ein ausführbares Minimalbeispiel angeben könntest. Am besten mit einem einfachen übersichtlichen Gleichungssystem $Au=b$ mit bekannter Lösung $u$ ohne Bezug auf das PDE-Problem.
Mir ist aufgefallen, dass Du die Abbruchbedingung (Zeile 5-7 in Algorithm 1) nicht implementiert hast.
Sollte außerdem $rohold$ nicht in der Schleife jeweils neu berechnet werden?
Weiter habe ich es mir nicht angeschaut.
Gruß
majoka
|
Profil
|
elliment hat die Antworten auf ihre/seine Frage gesehen. |
|
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]
|