Jacobi-metoden

TMA4400 Matematikk 1: Kalkulus og lineær algebra

Dato: 22. september 2025

Eksempel

La oss illustrere metoden ved å se på følgende lineære ligningssystem:

\[\begin{split} \left[\begin{array}{cc} 3 & 1 \\ 1 & 2 \end{array}\right]\cdot \left[\begin{array}{c} x \\ y \end{array}\right]=\left[\begin{array}{c} 5 \\ 5 \end{array}\right]. \end{split}\]

Vi ser at vi kan skrive

\[\begin{align*} x&=\frac{5-y}{3},\\ y&=\frac{5-x}{2}. \end{align*}\]

Med f. eks. startverdien \((x_0,y_0)=(0,0)\), gir dette oss følgende iterasjon:

\[\begin{align*} \left[\begin{array}{c} x_0 \\ y_0 \end{array}\right]&= \left[\begin{array}{c} 0 \\ 0 \end{array}\right]\\ \left[\begin{array}{c} x_1 \\ y_1 \end{array}\right]&= \left[\begin{array}{c} \frac{5-y_0}{3} \\ \frac{5-x_0}{2} \end{array}\right]=\left[\begin{array}{c} \frac{5}{3} \\ \frac{5}{2} \end{array}\right]\\ \left[\begin{array}{c} x_2 \\ y_2 \end{array}\right]&= \left[\begin{array}{c} \frac{5-y_1}{3} \\ \frac{5-x_1}{2} \end{array}\right]=\left[\begin{array}{c} \frac{5}{6} \\ \frac{5}{3} \end{array}\right]\\ &\vdots \end{align*}\]

Løsningen er \((x,y)=(1,2)\).

Jacobi-metoden

La \(A\) være ei \(n\times n\)-matrise og \(\mathbf{x},\mathbf{b}\in \mathbb{R}^n\) være vektorer. Vi skal nå formalisere algoritma vi så i eksempelet over. Vi skriver videre

\[ A=L+D+U, \]

hvor \(L\) er ei nedre triangulær matrise, \(D\) er ei diagonalmatrise og \(U\) ei øvre triangulær matrise. Legg merke til at notasjonen ikke er den samme som i \(LU\)-faktoriseringa, da var nemlig \(A=LU\).

Vi skal nå vise at det lineære ligningssystemet \(A\mathbf{x}=\mathbf{b}\) kan skrives på formen til en fikspunktiterasjon:

\[\begin{align*} A\mathbf{x}=\mathbf{b}\qquad&\Longleftrightarrow\qquad (D+L+U)\mathbf{x}=\mathbf{b},\\ D\mathbf{x}=\mathbf{b}-(L+U)\mathbf{x}\qquad&\Longleftrightarrow\qquad \mathbf{x}=D^{-1}(\mathbf{b}-(L+U)\mathbf{x}). \end{align*}\]

Altså har vi at \(\mathbf{x}=\mathbf{F}(\mathbf{x})\) der \(\mathbf{F}(\mathbf{x})=D^{-1}(\mathbf{b}-(L+U)\mathbf{x})\). Legg merke til at \(\mathbf{F}\) er en vektorfunksjon. Disse skal vi ikke analysere videre i dette emnet.

Jakobi-metoden kan nå beskrives som følgen \(\{\mathbf{a}_k\}_{k=0}^{\infty}\subseteq \mathbb{R}^n\) gitt rekursivt som

\[ \mathbf{a}_{k+1}=D^{-1}(\mathbf{b}-(L+U)\mathbf{a}_k),\qquad\text{med startverdien $\mathbf{a}_0$.} \]

Eksempel

Vi kan skrive om

\[\begin{split} \left[\begin{array}{cc} 3 & 1 \\ 1 & 2 \end{array}\right]\cdot \left[\begin{array}{c} x \\ y \end{array}\right]=\left[\begin{array}{c} 5 \\ 5 \end{array}\right], \end{split}\]

i eksempelet over, til en fikspunktiterasjon med \(\mathbf{a}_k=(x_k,y_k)\):

\[\begin{split} \left[\begin{array}{c} x_{k+1} \\ y_{k+1} \end{array}\right]=\mathbf{a}_{k+1}=D^{-1}(\mathbf{b}-(L+U)\mathbf{a}_k)=\left[\begin{array}{cc} 1/3 & 0 \\ 0 & 1/2 \end{array}\right]\left(\left[\begin{array}{c} 5 \\ 5 \end{array}\right]-\left(\left[\begin{array}{cc} 0 & 0 \\ 1 & 0 \end{array}\right]+\left[\begin{array}{cc} 0 & 1 \\ 0 & 0 \end{array}\right]\right)\cdot\left[\begin{array}{c} x_k \\ y_k \end{array}\right]\right)=\left[\begin{array}{c} \frac{5-y_k}{3} \\ \frac{5-x_k}{2} \end{array}\right], \end{split}\]

med startverdien \(\mathbf{a}_0=(0,0)\).

Vi implementerer nå metoden i python.

import numpy as np

# generisk iterativ fikspunktmetode
def jacobi(A, b, a, toleranse = 1.0E-10, K = 500):
    n = a.size               # antall elementer i a_0
    D = np.diag(A)*np.eye(n) # lager diagonalmatrisa
    Dinv = np.linalg.inv(D)  # lager inversen til diagonalmatrisa, dette er ikke dyrt siden den er diagonal
    L = np.tril(A,-1)        # lager L-matrisa
    U = np.triu(A,1)         # lager U-matrisa
    k = 0
    
    r = b-A@a # her bruker vi en annen notasjon for matrisemultiplikasjon når vi regner ut r
    feil = np.linalg.norm(r,ord=np.inf) # vi måler forskjellen mellom Ax og b i sup-normen
    
    while feil >= toleranse and k < K: # vi gjør dette mens ||Ax-b||<toleranse og k<K
        a = Dinv@(b-(L+U)@a) # ny startverdi
        k = k+1    # ny k
        r = b-A@a  # ny r
        feil = np.linalg.norm(r,ord=np.inf) # ny err
    if k >= K:
        print("Fant ikke fikspunkt for gitt startverdi. Returnerer nåværende verdier.")
        return a, feil, k
    else: 
        return a, feil, k

# definerer matrisa, høyresida og startvektoren vi skal bruke
A0 = np.array([[3.,1.],[1.,2.]])
b0 = np.array([[5.],[5.]]) # b som kolonnevektor 
a0 = np.array([[0.],[0.]]) # a0 som kolonnevektor

# finner fikspunkt med gitt startverdi
a, feil, k = jacobi(A0, b0, a0)

print(f"Løsning: {a}")
print(f"Feil: {feil}")
print(f"Iterasjoner: {k}")

Løsning: [[1.]
 [2.]]
Feil: 6.38049613144176e-11
Iterasjoner: 28

Eksempel

Vi skal nå se på det lineære ligningssystemet

\[\begin{split} \left[\begin{array}{cccc} 1 & 2 & 2 & 3 \\ -1 & 4 & 2 & 7 \\ 3 & 1 & 6 & 0 \\ 1 & 0 & 3 & 4 \end{array}\right]\cdot \left[\begin{array}{c} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}\right]=\left[\begin{array}{c} 0 \\ 1 \\ -1 \\ 2 \end{array}\right]. \end{split}\]

Dette systemet har løsningen \((x_1,x_2,x_3,x_4)=(-0.1979, -0.7187, 0.0521, 0.5104)\), men Jacobi-metoden finner den ikke:

# definerer matrisa, høyresida og startvektoren vi skal bruke
A1 = np.array([[1.,2.,2.,3.],[-1.,4.,2.,7.],[3.,1.,6.,0],[1.,0.,3.,4.]])
b1 = np.array([[0.],[1.],[-1.],[2.]]) # b som kolonnevektor 
a1 = np.array([[0.],[0.],[0.],[0.]]) # a0 som kolonnevektor

# finner fikspunkt med gitt startverdi
a, feil, k = jacobi(A1, b1, a1)

print(f"Løsning: {a}")
print(f"Feil: {feil}")
print(f"Iterasjoner: {k}")

Fant ikke fikspunkt for gitt startverdi. Returnerer nåværende verdier.
Løsning: [[-7.23872269e+125]
 [-1.45973794e+125]
 [-2.15594313e+125]
 [-1.91258198e+125]]
Feil: 3.6111564766014174e+126
Iterasjoner: 500

Tilstrekkelig betingelse for konvergens

Vi sier at ei \(n\times n\)-matrise \(A\) med elementer \(\{a_{ij}\}_{i,j}\) er strengt diagonaldominant dersom

\[ |a_{ii}|>\sum_{j\neq i}|a_{ij}|,\qquad\text{for alle $i=1,\ldots,n$.} \]

---

Teorem. La \(A\) være ei strengt diagonaldominant \(n\times n\)-matrise. Da konvergerer Jacobi-metoden for ligningssystemet \(A\mathbf{x}=\mathbf{b}\) mot en entydig løsning \(\mathbf{x}\in\mathbb{R}^n\) av systemet for alle gitte høyresider \(\mathbf{b}\in\mathbb{R}^n\) og alle startverdier \(\mathbf{a}_0\in\mathbb{R}^n\).

---

Eksempel

Vi skal nå se på det lineære ligningssystemet

\[\begin{split} \left[\begin{array}{ccc} 10 & -1 & 2 \\ -1 & 11 & -1 \\ 2 & -1 & 10 \end{array}\right]\cdot \left[\begin{array}{c} x_1 \\ x_2 \\ x_3 \end{array}\right]=\left[\begin{array}{c} 6 \\ 22 \\ -10 \end{array}\right]. \end{split}\]

La oss prøve å løse dette med Jacobi-metoden:

# definerer matrisa, høyresida og startvektoren vi skal bruke
A2 = np.array([[10.,-1.,2.],[-1.,11.,-1.],[2.,-1.,10.]])
b2 = np.array([[6.],[22.],[-10.]]) # b som kolonnevektor 
a2 = np.array([[0.],[0.],[0.]]) # a0 som kolonnevektor

# finner fikspunkt med gitt startverdi
a, feil, k = jacobi(A2, b2, a2)

print(f"Løsning: {a}")
print(f"Feil: {feil}")
print(f"Iterasjoner: {k}")

Løsning: [[ 1.]
 [ 2.]
 [-1.]]
Feil: 2.7418067816142866e-11
Iterasjoner: 20

Eksempel

Vi skal nå gå tilbake til det første lineære ligningssystemet vi så på

\[\begin{split} \left[\begin{array}{cc} 3 & 1 \\ 1 & 2 \end{array}\right]\cdot \left[\begin{array}{c} x \\ y \end{array}\right]=\left[\begin{array}{c} 5 \\ 5 \end{array}\right]. \end{split}\]

Dette var strengt diagonaldominant. Hvis vi bytter om på første og andre rad, får vi

\[\begin{split} \left[\begin{array}{cc} 1 & 2 \\ 3 & 1 \end{array}\right]\cdot \left[\begin{array}{c} x \\ y \end{array}\right]=\left[\begin{array}{c} 5 \\ 5 \end{array}\right]. \end{split}\]

Dette systemet er imidlertid ikke strengt diagonaldominant.

# definerer matrisa, høyresida og startvektoren vi skal bruke
A3 = np.array([[1.,2.],[3.,1]])
b3 = np.array([[5.],[5.]]) # b som kolonnevektor 
a3 = np.array([[0.],[0.]]) # a0 som kolonnevektor

# finner fikspunkt med gitt startverdi
a, feil, k = jacobi(A3, b3, a3)

print(f"Løsning: {a}")
print(f"Feil: {feil}")
print(f"Iterasjoner: {k}")

Fant ikke fikspunkt for gitt startverdi. Returnerer nåværende verdier.
Løsning: [[-3.44994837e+194]
 [-6.89989674e+194]]
Feil: 1.7249741857848466e+195
Iterasjoner: 500