Matroids Matheplanet Forum Index
Moderiert von Bilbo
Informatik » Angewandte Informatik » cg solver in c++
Autor
Universität/Hochschule cg solver in c++
elliment
Neu Letzter Besuch: vor mehr als 3 Monaten
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 Letzter Besuch: in der letzten Woche
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 Letzter Besuch: in der letzten Woche
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 Letzter Besuch: in der letzten Woche
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.

Wechsel in ein anderes Forum:
 Suchen    
 
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]