Iterative Lösung nichtlinearer gleichungen

Gast #7316491

Ich habe eine Gleichung x = g(x) die nicht algebraisch für x lösbar ist. Wie kann man am besten in c++ solche Gleichungen iterativ lösen? Ich hatte mich mal mit der dlib beschäftigt und im Prinzip muss man doch x - g(x) = 0 lösen, oder? Wie würde denn die objektive Funktion dafür aussehen?

Persönliche Seite #7316961

Erste Methode ist Newton-Raphson Iteration für ƒ(x) = g(x) - x:

  1. Startwert für Iteration raten oder ausarbeiten

  2. Iterieren bis zur gewünschten Genauigkeit:

$$x\mapsto x-\frac{f(x)}{f'(x)} = x-\frac{g(x)-x}{g'(x)-1} = \frac{xg'(x)-g(x)}{g'(x)-1}$$

  • Funktioniert gut (quadratische Konvergenz) für einfache Nullstellen falls man die Iteration mit guter[tm] Näherung startet.

  • Schlechtere Konvergenz (linear) bei mehrfacher Nullstelle.

  • Es gibt Startwerte, für welche die Iteration nicht konvergiert, und diese können eine Menge von vollem Maß bilden.

Gast #7317401

Maxe schrieb:

Carol Z. schrieb:

g(x) = a * 1 / (1 - (1 + sqrt(1 / x))^-4) also x - g(x) = 0

Das "also" macht hier keinen Sinn. Bist du sicher dass g(x) = x ist?

Ja da bin ich mir sicher. Wie gesagt das ist das Optimierungsproblem für x, das es zu lösen gilt.

Moderator #7317443

Yalu X. schrieb:

Für a≠0 hat die Funktion keine reellen Nullstellen.

Es geht nicht um die Nullstellen von g(x) sondern von x-g(x), da wird es schon Nullstelle(n) geben.

Upps, das hatte ich übersehen. Ja du hast recht, dann gibt es natürlich (mindestens) eine reelle Lösung.

Da die Funktion fest vorgegeben ist, kann man das mit Newton lösen, allerdings ist dann die Festlegung des Startwerts nicht ganz trivial. Dasselbe Problem hat man auch mit der Regula falsi.

Man kann die Gleichung aber auch in eine algebraische Gleichung in √x umformen. Dafür gibt es Algorithmen, die ohne die Vorgabe eines Startwerts sämtliche Lösungen liefern, bspw. den PolynomSolver in eigen3/unsupported. Das sieht dann etwa so aus:

solve.cpp:

1
#include <iostream>
2
#include <cstdlib>
3
#include <eigen3/unsupported/Eigen/Polynomials>
4

5
int main(int argc, char *argv[]) {
6
  const int degree = 5;
7
  const double a = strtod(argv[1], NULL);
8
  Eigen::Vector<double, degree+1> coeffs(a, 4*a, 6*a-1, 4*a-4, a-6, -4);
9
  Eigen::PolynomialSolver<double, degree> solver;
10
  solver.compute(coeffs);
11
  const auto roots = solver.roots();
12
  for (auto c: roots) {
13
    const auto realpart = c.real();
14
    if (c.imag() == 0 && realpart > 0) {
15
      const auto solution = realpart * realpart;
16
      std::cout << solution << '\n';
17
    }
18
  }
19
}

Anwendungsbeispiel:

1
$ solve 4.2  # a = 4.2
2
5.54664

Gegenprobe:

1
4.2 * 1 / (1 - (1 + sqrt(1 / 5.54664))**-4) = 5.54664

Neben der Eigen-Bibliothek gibt es noch jede Menge anderer Bibliotheken, die eine ähnliche Funktion bereitstellen. Am besten suchst du einfach nach Polynomial Solver.

Antwort schreiben

Eine eigene Frage? Dieser Thread ist über 3 Jahre alt. Stell sie in einem eigenen Thread – dort lesen mehr Leute mit, und sie wird öfter und schneller beantwortet.

Neuen Thread starten

Eine Lösung, Ergänzung oder Rückmeldung zu diesem Thread? Dann antworte hier.

Formatierung
  • [c]C-Code[/c]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] – Link zu einem Artikel
  • Verweis auf einen anderen Beitrag: Rechtsklick auf dessen Beitrags-ID (#8078628), „Adresse kopieren“, in den Text einfügen

Mehr zur Formatierung

Längeren Sourcecode als Anhang, nicht in den Text. JPG nur für Fotos und Scans; Zeichnungen und Screenshots als PNG oder GIF. Siehe Bildformate.

Schon registriert? Anmelden

Wird nicht angezeigt. Wir schicken dir einen Link; sobald du damit bestätigst, erscheint dein Beitrag.

Mit dem Abschicken legst du einen kostenlosen Account an und erkennst die Nutzungsbedingungen und die Datenschutzerklärung an.