シェアする

  • このエントリーをはてなブックマークに追加

How to solve simultaneous linear equation with a tridiagonal matrix

シェアする

  • このエントリーをはてなブックマークに追加
猫の絵

This time I will show you how to solve simultaneous linear equations with a tridiagonal matrix.

The content is as follows.

I also have other articles here about simulation.

Tridiagonal matrix

In linear algebra, a tridiagonal matrix is a band matrix that has nonzero elements on the main diagonal, the first diagonal below this, and the first diagonal above the main diagonal only. (Wikipedia)

\(\begin{eqnarray} A = \left( \begin{array}{ccccc} b_{1} & c_{1} & & & \\ a_{2} & b_{2} & c_{2} & & \\ & a_{3} & b_{3} & c_{3} & \\ & & \ddots & \ddots & \\ & & & a_{n} & b_{n} \\ \end{array} \right) \end{eqnarray}\)

All elements in the empty space are zero.
A simultaneous linear equation with a tridiagonal matrix as the coefficient matrix can be solved by \(O\left(n\right)\) (\(n\) is the matrix size).
I will introduce it here.

The way to solve

We consider solving a simultaneous linear equation

\(\begin{eqnarray} A{\bf x} = {\bf d}, \end{eqnarray}\)

here, \({\bf x},\,{\bf d}\) are both \(n\) dimensional vectors.

Jumping right to the conclusion, with \(P_{i},\,Q_{i}\) that satisfy

\(\begin{eqnarray} P_{1} &=& -\frac{c_{1}}{a_{1}} \\ Q_{1} &=& \frac{d_{1}}{a_{1}} \\ P_{i} &=& -\frac{c_{i}}{b_{i}+a_{i}P_{i-1}} \\ Q_{i} &=& \frac{d_{i}+a_{i}Q_{i-1}}{b_{i}+a_{i}P_{i-1}} \end{eqnarray}\)

it can be solved as

\(\begin{eqnarray} x_{n} &=& Q_{n} \\ x_{i} &=& P_{i}x_{i+1}+Q_{i} \end{eqnarray}\)

First, we look at the first line.
From

\(\begin{eqnarray} b_{1}x_{1}+c_{1}x_{2} = d_{1}, \end{eqnarray}\)

we get

\(\begin{eqnarray} x_{1} = -\frac{c_{1}}{b_{1}}x_{2}+\frac{d_{1}}{b_{1}} \end{eqnarray}\)

Now, for later use, let \( P_{1} = -c_{1}/b_{1},\,Q_{1} = d_{1}/b_{1}\), then we can write

\(\begin{eqnarray} x_{1} = P_{1}x_{2}+Q_{2} \end{eqnarray}\)

Next, the second line.

\(\begin{eqnarray} a_{2}x_{1}+b_{2}x_{2}+c_{2}x_{3} = d_{2} \end{eqnarray}\)

Substituting the above \(x_{1}\) and organizing,

\(\begin{eqnarray} x_{2} = P_{2}x_{3}+Q_{2}, \end{eqnarray}\)

here,

\(\begin{eqnarray} P_{2} &=& -\frac{c_{2}}{a_{2}P_{1}+b_{2}} \\ Q_{2} &=& \frac{d_{2}-a_{2}Q_{1}}{a_{2}P_{1}+b_{2}} \end{eqnarray}\)

Repeating the same thing, up to the \(\left(n-1\right)\)-th line,

\(\begin{eqnarray} x_{i} &=& P_{i}x_{i+1}+Q_{i} \\ P_{i} &=& -\frac{c_{i}}{a_{i}P_{i-1}+b_{i}} \\ Q_{i} &=& \frac{d_{i}-a_{i}Q_{i-1}}{a_{i}P_{i-1}+b_{i}} \end{eqnarray}\)

Lastly, the \(n\)-th line is

\(\begin{eqnarray} a_{n}x_{n-1}+b_{n}x_{n} = d_{n}, \end{eqnarray}\)

but the equation of the \(\left(n-1\right)\)-th line

\(\begin{eqnarray} x_{n-1} &=& P_{n-1}x_{n}+Q_{n-1} \\ \end{eqnarray}\)

gives us the following result,

\(\begin{eqnarray} x_{n} = Q_{n} \\ x_{n-1} = P_{n-1}x_{n}+Q_{n-1} \end{eqnarray}\)

After all, we can get the answer \(x_{i}\) with \(O\left(n\right)\) calculations by finding \(P_{i},\,Q_{i}\) in the previous calculation.

Implementation example

void trimat_sol(double **A,double *a,double *u,int n){
  A[0][1] = -A[0][1];
  for(int i = 1;i<n-1;i++){
    A[i][i-1] = -A[i][i-1];
    A[i][i+1] = -A[i][i+1];
  }
  A[n-1][n-2] = -A[n-1][n-2];

  double *P,*Q;
  if((P = (double *)malloc(sizeof(double)*n))==NULL){
    printf("*** Allocation ERROR P ***\n");
    exit(1);
  }
  if((Q = (double *)malloc(sizeof(double)*n))==NULL){
    printf("*** Allocation ERROR Q ***\n");
    exit(1);
  }
  
  double den;
  double num;
  P[0] = A[0][1]/A[0][0];
  Q[0] = a[0]/A[0][0];
  for(int i = 1;i<n;i++){
    den = A[i][i]-A[i][i-1]*P[i-1];
    P[i] = A[i][i+1]/den;
    num = a[i]+A[i][i-1]*Q[i-1];
    Q[i] = num/den;
  }

  u[n-1] = Q[n-1];
  for(int i = n-2;i>=0;i--){
    u[i] = P[i]*u[i+1]+Q[i];
  }
  free(P);
  free(Q);
}

An exercise

Given matrix \(A\), solve the equation

\(\begin{eqnarray} A{\bf x} = {\bf d} \end{eqnarray}\)

with \(A\) and \(d\)

\(\begin{eqnarray} A &=& \left( \begin{array}{cccc} 2 & 1 & 0 & 0 \\ 1 & 2 & 1 & 0 \\ 0 & 1 & 2 & 1 \\ 0 & 0 & 1 & 2 \\ \end{array} \right) \\ {\bf d} &=& \left( \begin{array}{c} 3 \\ 6 \\ 9 \\ 10 \\ \end{array} \right) \end{eqnarray}\)

The answer;

\(\begin{eqnarray} {\bf x} = \left( \begin{array}{c} 0.4 \\ 2.2 \\ 1.2 \\ 4.4 \\ \end{array} \right) \end{eqnarray}\)

Summary

Now, that’s all what I wanted to write.

This time, I introduced a way to solve a simultaneous linear equation with a tridiagonal coefficient matrix.

Thank you for reading and please spread this blog if you like.

If you have any comment, please let me know from the e-mail address below or the CONTACT on the menu bar.
tsunetthi(at)gmail.com
Please change (at) to @.

Or you can also contact me via twitter (@warotan3)

シェアする

  • このエントリーをはてなブックマークに追加

フォローする