Mathematische Begründung

Es wird eine Faktorisierung ${\displaystyle f(x)=a(x)b(x)}$ mit einem quadratischen Faktor a(x) und einem verbleibenden Faktor b(x) gesucht. Ist jedoch a(x) nur ein näherungsweiser Faktor, so hinterlässt die Division mit Rest von f(x) durch a(x) mit Ergebnis b(x) einen Rest r(x),

${\displaystyle f(x)=a(x)b(x)+r(x)}$.

Da a(x) quadratisch ist, ist r(x) linear oder konstant. Es wird nun ein lineares Polynom ${\displaystyle \Delta a(x)}$ gesucht, sodass ${\displaystyle a(x)+\Delta a(x)}$ eine bessere Näherung für einen Faktor von f(x) ist. Ist das verbesserte Polynom sogar ein exakter Faktor, so gibt es ein Polynom ${\displaystyle \Delta b(x)}$ vom Grad n-3 mit

${\displaystyle f(x)=\left(a(x)+\Delta a(x)\right)\cdot \left(b(x)+\Delta b(x)\right)}$

Im Newton-Verfahren werden nur Terme erster Ordnung in den Koordinaten des Änderungsvektors betrachtet, alle höhergradigen Ausdrücke werden vernachlässigt. In erster Ordnung ergibt das unter Kombination beider Gleichungen

${\displaystyle r(x)=f(x)-a(x)b(x)=a(x)\cdot \Delta b(x)+\Delta a(x)\cdot b(x)}$.

Diese Gleichung kann in ein lineares Gleichungssystem für die Koeffizienten von ${\displaystyle \Delta a(x)}$ und ${\displaystyle \Delta b(x)}$ aufgelöst werden. Die Gleichung kann jedoch mit algebraischen Mitteln weiter vereinfacht werden, so dass das schließlich zu lösende lineare System zwei Variable und genauso viele Gleichungen hat.

Da ${\displaystyle \deg \Delta a(x)<\deg a(x)}$ gilt, ist das Polynom ${\displaystyle \Delta a(x)}$ schon selbst der Vertreter kleinsten Grades in seiner Restklasse modulo a(x). Da ${\displaystyle a(x)\Delta b(x)}$ in der Restklasse der Null liegt, muss

${\displaystyle r(x)\equiv \Delta a(x)\,b(x)\mod a(x)}$

gelten. Nun kann auch das Polynom b(x) modulo a(x) reduziert werden, nach einer weiteren Division mit Rest ergeben sich ein Quotient q(x) und ein linearer Rest p(x) mit ${\displaystyle b(x)=q(x)a(x)+p(x)}$. Es ergibt sich

${\displaystyle r(x)\equiv \Delta a(x)p(x)\mod a(x)}$ bzw. ${\displaystyle r(x)=\Delta a(x)p(x)-p_{1}\Delta a_{1}\cdot a(x)}$.

Als Gleichungssystem ausgeschrieben bedeutet dies

${\displaystyle {\begin{pmatrix}r_{0}\\r_{1}\end{pmatrix}}={\begin{pmatrix}p_{0}&-p_{1}a_{0}\\p_{1}&p_{0}-p_{1}a_{1}\end{pmatrix}}\cdot {\begin{pmatrix}\Delta a_{0}\\\Delta a_{1}\end{pmatrix}}}$.

Die Determinante der Systemmatrix ist ${\displaystyle D=p_{0}^{2}+p_{1}(a_{0}p_{1}-a_{1}p_{0})}$. Ist diese von Null verschieden, so ergibt sich die Lösung des Systems und damit die Änderung des quadratischen Faktors g(x) als

${\displaystyle {\begin{pmatrix}\Delta a_{0}\\\Delta a_{1}\end{pmatrix}}={\frac {1}{D}}{\begin{pmatrix}p_{0}-p_{1}a_{1}&p_{1}a_{0}\\-p_{1}&p_{0}\end{pmatrix}}\cdot {\begin{pmatrix}r_{0}\\r_{1}\end{pmatrix}}={\begin{pmatrix}p_{0}r_{0}+p_{1}(a_{0}r_{1}-a_{1}r_{0})\\p_{0}r_{1}-p_{1}r_{0}\end{pmatrix}}}$.

Für den nächsten Schritt wird a(x) durch ${\displaystyle a(x)+\Delta a(x)}$ ersetzt und von vorn begonnen. Man kann abbrechen, wenn die Koeffizienten des Rests r(x) sämtlich eine vorher gesetzte Schranke unterschreiten.

Grundzüge

Das Verfahren beruht auf folgenden Schritten:

1. Aus den Startwerten der Iteration für zwei Nullstellen wird ein quadratisches Polynom konstruiert.
2. Das Polynom ${\displaystyle n}$-ten Grades wird durch dieses quadratische Polynom dividiert
3. Das bei der Division entstehende Polynom ${\displaystyle (n-2)}$-ten Grades wird erneut durch das quadratische Polynom dividiert
4. Bei der Division entstehen Divisionsreste
5. Aus den Divisionsresten wird ein verbessertes quadratisches Polynom bestimmt
6. Wenn bei der Polynomdivision kein Rest mehr bleibt, sind die Nullstellen des quadratischen Polynoms auch Nullstellen des Polynoms n-ten Grades.

Iteration

• Durch Polynomdivision mit dem quadratischen Polynom

${\displaystyle a(x)=x^{2}+a_{1}\cdot x+a_{0}}$

erhält man ein Polynom ${\displaystyle (n-2)}$-ten Grades:

${\displaystyle b(x)=b_{n-2}\cdot x^{n-2}+b_{n-3}\cdot x^{n-3}+\dots +b_{1}\cdot x+b_{0}}$

Die Koeffizienten des zugehörigen Divisionsrests ${\displaystyle r(x)=r_{0}+r_{1}x=p-ab}$ können mit der folgenden Rekursionsformel zur Bestimmung von ${\displaystyle \,b_{i}}$ berechnet werden:

${\displaystyle b_{n}=b_{n-1}=0}$;

${\displaystyle b_{j}=p_{j+2}-a_{1}\cdot b_{j+1}-a_{0}\cdot b_{j+2};\qquad \qquad j=n-2,n-3,\dots ,0,-1,-2}$;

${\displaystyle r_{1}=b_{-1}}$ und

${\displaystyle r_{0}=b_{-2}+a_{1}b_{-1}}$;

• Erneute Division von b(x) durch a(x) ergibt ein neues Polynom ${\displaystyle q(x)}$ und dessen Divisionsrest ${\displaystyle p(x)=p_{1}x+p_{0}\!}$, wobei eine analoge Rekursionsformel zum Einsatz kommt, es gelten wieder ${\displaystyle p_{1}=q_{-1}}$ und ${\displaystyle p_{0}=q_{-2}+a_{1}q_{-1}}$;

• Mit den Divisionsresten ${\displaystyle r_{0},\ r_{1},\ p_{0},\ p_{1}}$ verbessert man die Koeffizienten des quadratischen Polynoms a(x):

${\displaystyle M=-a_{0}\cdot q_{-1}-a_{1}\cdot q_{-2};\qquad D=q_{-2}^{2}-M\cdot q_{-1}}$

${\displaystyle \Delta a_{1}={\frac {q_{-1}\cdot b_{-2}-q_{-2}\cdot b_{-1}}{D}};\qquad \Delta a_{0}={\frac {M\cdot b_{-1}-q_{-2}\cdot b_{-2}}{D}}}$

• Die verbesserten Startwerte ergeben sich aus

${\displaystyle a_{1}^{neu}=a_{1}^{alt}-\Delta a_{1}\,}$

${\displaystyle a_{0}^{neu}=a_{0}^{alt}-\Delta a_{0}\,}$

Anmerkungen zum Ablauf

• Die Bestimmung von ${\displaystyle \Delta a_{0},\ \Delta a_{1}}$ kann als Lösung des folgenden Gleichungssystems aufgefasst werden:

${\displaystyle {\begin{pmatrix}q_{-2}&q_{-1}\\M&q_{-2}\end{pmatrix}}\cdot {\begin{pmatrix}\Delta a_{1}\\\Delta a_{0}\end{pmatrix}}={\begin{pmatrix}-b_{-1}\\-b_{-2}\end{pmatrix}}}$

• Das Verfahren lässt sich numerisch relativ einfach, schnell und speichergünstig umsetzen, wenn man nicht erst ${\displaystyle \,b(x)}$ berechnet und speichert, sondern im Schritt j die ${\displaystyle \,b_{j},b_{j+1},b_{j+2}}$ sofort dazu verwendet, auch den Koeffizienten ${\displaystyle \,q_{j}=b_{j+2}-a_{1}q_{j+1}-a_{0}q_{j+2}}$ zu berechnen.

• Als Startwerte für die Iteration kann man beispielsweise aus den führenden drei Koeffizienten das quadratische Polynom ${\displaystyle a(x)={\frac {1}{f_{n}}}(f_{n}x^{2}+f_{n-1}+f_{n-2}}$ bilden, d. h. ${\displaystyle a_{1}={\frac {f_{n-1}}{f_{n}}};\quad a_{0}={\frac {f_{n-2}}{f_{n}}}}$ wählen. Wenn Näherungen für zwei Nullstellen ${\displaystyle {\tilde {x}}_{1},\ {\tilde {x}}_{2}}$ vorliegen, kann man auch alternativ nach den Vietaschen Formeln ${\displaystyle a_{1}=-{\tilde {x}}_{1}-{\tilde {x}}_{2};\quad a_{0}={\tilde {x}}_{1}\cdot {\tilde {x}}_{2}}$ wählen.

• Die Iteration kann abgebrochen werden, wenn ${\displaystyle \,r_{0}=r_{1}=0}$ bzw. ${\displaystyle \,b_{-1}=b_{-2}=0}$ in der gewünschten Ergebnisgenauigkeit vorliegen, dann wurden passende Nullstellen gefunden und ${\displaystyle \,a(x)}$ ist ein Teiler von ${\displaystyle \,f(x)}$. In diesem Fall bestimmt man alle Koeffizienten von ${\displaystyle \,b(x)}$ und dann mit dem neuen Polynom ${\displaystyle \,f(x)=b(x)}$ reduzierten Grades die nächsten Nullstellen zu suchen. Denn wenn man von ${\displaystyle \,f(x)}$ die Linearfaktoren zu den Nullstellen abdividiert, erhält man ${\displaystyle \,b(x)}$.

Algorithmus

Das Programmbeispiel in Pseudocode beschreibt einen einzelnen Iterationsschritt, bei dem die Koeffizienten a0 und a1 des quadratischen Polynoms durch zwei Korrekturen da0 und da1 verbessert werden. Das Feld (Array) f() enthält das Polynom, wobei in f(n) der Koeffizient von ${\displaystyle x^{n}}$ steht. Die Bezeichnungen der Variablen wurden sonst wie in den Formeln gewählt.

In der Schleife der doppelten Polynomdivision stehen die Programmvariablen b0, b1, b2 im Schritt j für die Koeffizienten ${\displaystyle \,b_{j-2},b_{j-1},b_{j}}$ und entsprechend q0, q1, q2 für ${\displaystyle \,q_{j-2},q_{j-1},q_{j}}$. Beim Übergang von Schritt j zu Schritt j−1 müssen die Programmvariablen einen Koeffizientenindex weiter geschoben werden, b2←b1←b0 und q2←q1←q0.

    j=n
    b0 = f(n)
    b1 = 0
    q0 = 0
    q1 = 0

    For j = n - 1 To 0 Step -1
        b2 = b1
        b1 = b0
        b0 = f(j) - a1 * b1 - a0 * b2

        q2 = q1
        q1 = q0
        q0 = b2   - a1 * q1 - a0 * q2
    Next j

    M  =  -a0 * q1 - a1  * q0
    D  =  q0 * q0 - M  * q1
    da1 = (b0 * q1 - b1 * q0) / D
    da0 = (b1 * M  - b0 * q0) / D

    a1  = a1 - da1
    a0  = a0 - da0

Im Ergebnis des Iterationsschritts enthalten die Programmvariablen a0 und a1 verbesserte Werte für die Koeffizienten von ${\displaystyle \,a(x)=x^{2}+a_{1}x+a_{0}}$