# How to obtain a convergent solution iteratively for a linear system of equations? [closed]

I am working on a problem that requires an iterative procedure to solve a linear system of equations, the system of equations in matrix form is:

$$\underbrace{\begin{bmatrix} r_{11} & r_{12} & r_{13} & \cdots & r_{1j} \\ r_{21} & r_{22} & r_{23} & \cdots & r_{2j} \\ r_{31} & r_{32} & r_{33} & \cdots & r_{3j} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ r_{i1} & r_{i2} & r_{i3} & \cdots & r_{ij} \end{bmatrix}}_\textit{R} \underbrace{\begin{bmatrix} a_{1} & 0 & 0 & \cdots & 0 \\ 0 & a_{2} & 0 & \cdots & 0 \\ 0 & 0 & a_{3} & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & a_{j} \end{bmatrix}}_\textit{A} + \underbrace{\begin{bmatrix} b_{1} & b_{2} & b_{3} & \cdots & b_{j} \\ b_{1} & b_{2} & b_{3} & \cdots & b_{j} \\ b_{1} & b_{2} & b_{3} & \cdots & b_{j} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ b_{1} & b_{2} & b_{3} & \cdots & b_{j} \end{bmatrix}}_\textit{B} = \underbrace{\begin{bmatrix} c_{1} & c_{1} & c_{1} & \cdots & c_{1} \\ c_{2} & c_{2} & c_{2} & \cdots & c_{2} \\ c_{3} & c_{3} & c_{3} & \cdots & c_{3} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ c_{i} & c_{i} & c_{i} & \cdots & c_{i} \end{bmatrix}}_\textit{C}\\ R_{i, j}A_{j, j} + B_{i, j} = C_{i, j}$$ Now matrix R is fully known (input), matrices A, B, C are unknown. I am working on an iterative procedure in which I provide initial guesses for A and B, then calculate C. The iteration is to be carried out till I get a converged value for C (output I require). The code is being developed in Python. So far no luck in choosing the initial guess too.

EXAMPLE

To illustrate the problem I have given an example for a $2 \times 2$ matrix:

$$\begin{bmatrix} r_{11} & r_{12} \\ r_{21} & r_{22} \\ \end{bmatrix} \begin{bmatrix} a_{1} & 0 \\ 0 & a_{2} \\ \end{bmatrix} + \begin{bmatrix} b_{1} & b_{2} \\ b_{1} & b_{2} \\ \end{bmatrix} = \begin{bmatrix} c_{1} & c_{1} \\ c_{2} & c_{2} \\ \end{bmatrix}\\$$ this gives us:

$$a_{1}r_{11} + b_{1} = a_{2}r_{12} + b_{2} = c_{1}\\ a_{1}r_{21} + b_{1} = a_{2}r_{22} + b_{2} = c_{2}$$

From these equations I get:

$$a_{2} = a_{1}\frac{(r_{11} - r_{21})}{(r_{12} - r_{22})}; \quad b_{2}=b_{1} + a_{1}\frac{(r_{11} - r_{21})(r_{11} - r_{12})}{(r_{12} - r_{22})}$$

This helps in the iterative method, but I would like to generalize it for a bigger matrix.

UPDATE

I have added two reproducible Python codes that should be helpful for a $2 \times 2$ matrix:

### Code 1

import numpy as np
R = np.matrix([[2.5, 2.9], [2.3, 2.7]])
m = R.shape[0]
n = R.shape[1]
A = np.matrix(np.zeros(shape=(n,n)))
B = np.matrix(np.zeros(shape=(m,n)))
A[0,0] = 1
B[0,0] = 1
k = 0
for i in range (1, 10000):
A[0,0] = A[0, 0] - 0.0000001
B[0,0] = B[0, 0] - 0.0000001
for j in range(1, n):
A[j, j] = A[(j-1), (j-1)] * ((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1]))
B[0, j] = B[0, 0] + A[0, 0] * (((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1])) * (R[0, 0] - R[0, 1]))
B[1, 0] = B[0, 0]
B[1, j] = B[0, 0] + A[0, 0] * (((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1])) * (R[0, 0] - R[0, 1]))
C = R * A + B
C_convergence = np.all(np.all(np.diff(C, axis = 1) == 0, axis = 1), axis = 0)
k = k + 1
print k
print C


### Code 2

import numpy as np
R = np.matrix([[2.5, 2.9], [2.3, 2.7]])
m = R.shape[0]
n = R.shape[1]
A = np.matrix(np.zeros(shape=(n,n)))
B = np.matrix(np.zeros(shape=(m,n)))
A[0,0] = 1
B[0,0] = 1
k = 0
C_convergence = False
k = 0
while C_convergence == False:
A[0,0] = A[0, 0] - 0.0000001
B[0,0] = B[0, 0] - 0.0000001
for j in range(1, n):
A[j, j] = A[(j-1), (j-1)] * ((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1]))
B[0, j] = B[0, 0] + A[0, 0] * (((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1])) * (R[0, 0] - R[0, 1]))
B[1, 0] = B[0, 0]
B[1, j] = B[0, 0] + A[0, 0] * (((R[0, 0] - R[1, 0]) / (R[0, 1] - R[1, 1])) * (R[0, 0] - R[0, 1]))
C = R * A + B
C_convergence = np.all(np.all(np.diff(C, axis = 1) == 0, axis = 1), axis = 0)
k = k + 1
print k
print C


## UPDATE 2

I have written the matrices in an alternative way which is of the form $AX = B$:

$$\begin{bmatrix} r_{11} & 0 & 0 & \cdots & 0 & 1 & 0 & 0 & \cdots & 0 & -1 & 0 & 0 & \cdots & 0 \\ 0 & r_{12} & 0 & \cdots & 0 & 0 & 1 & 0 & \cdots & 0 & -1 & 0 & 0 & \cdots & 0 \\ 0 & 0 & r_{13} & \cdots & 0 & 0 & 0 & 1 & \cdots & 0 & -1 & 0 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & r_{1j} & 0 & 0 & 0 & \cdots & 1 & -1 & 0 & 0 & \cdots & 0 \\ r_{21} & 0 & 0 & \cdots & 0 & 1 & 0 & 0 & \cdots & 0 & 0 & -1 & 0 & \cdots & 0 \\ 0 & r_{22} & 0 & \cdots & 0 & 0 & 1 & 0 & \cdots & 0 & 0 & -1 & 0 & \cdots & 0 \\ 0 & 0 & r_{23} & \cdots & 0 & 0 & 0 & 1 & \cdots & 0 & 0 & -1 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & r_{2j} & 0 & 0 & 0 & \cdots & 1 & 0 & -1 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots &\\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots &\\ r_{i1} & 0 & 0 & \cdots & 0 & 1 & 0 & 0 & \cdots & 0 & 0 & 0 & 0 & \cdots & -1 \\ 0 & r_{i2} & 0 & \cdots & 0 & 0 & 1 & 0 & \cdots & 0 & 0 & 0 & 0 & \cdots & -1 \\ 0 & 0 & r_{i3} & \cdots & 0 & 0 & 0 & 1 & \cdots & 0 & 0 & 0 & 0 & \cdots & -1 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & r_{ij} & 0 & 0 & 0 & \cdots & 1 & 0 & 0 & 0 & \cdots & -1 \\ \end{bmatrix} \begin{bmatrix} a_{1} \\ a_{2} \\ a_{3} \\ \vdots \\ a_{j} \\ b_{1} \\ b_{2} \\ b_{3} \\ \vdots \\ b_{j} \\ c_{1} \\ c_{2} \\ c_{3} \\ \vdots \\ c_{i} \\ \end{bmatrix} = 0$$

Now how do I solve this using an iterative method?

• Just a stab in the dark, You could probably write this as a convex optimization problem where you minimize some matrix norm of A,B,C subject to the relation RA+B=C. I know that python has some packages that do convex optimization. – Baby Dragon Feb 22 '15 at 6:00
• I'm voting to close this question as off-topic because it makes no sense mathematically. – David Ketcheson Apr 24 '15 at 8:03
• The problem is vastly underdetermined. As one trivial solution, take $A=0$ and $B=C$ where $C$ is arbitrary. – David Ketcheson Apr 24 '15 at 8:04
• I agree with David, there's no hope without further constraints on A,B,C. You should edit the question to include some information about the background of this problem (always a good idea!), maybe some inherent constraints can be found. Otherwise you need a completely different approach. In the meantime, I am voting to put this question on hold (and will vote to reopen once edited) because it has no useful answer in its current form. – Christian Clason Apr 24 '15 at 8:31

It looks like you have i total equations. If A, B, and C are all unknown then you have j + j + i number of unknowns. You cannot find a unique solution to this problem because you have more unknowns than number of equations. Even if you create an iterative method that converges to a solution, this will be one solution of an infinite number of solutions.

However, for people interested in iterative solutions to linear systems of equations they can look into the following:

Jacobi Method (http://en.wikipedia.org/wiki/Jacobi_method)

Gauss-Seidel (http://en.wikipedia.org/wiki/Gauss-Seidel_method)

For nonlinear algebraic equations you can look into the Newton-Raphson method.

• FYI there would be i x j total equations. – nxkryptor Feb 20 '15 at 18:31
• @nxkryptor My mistake, thanks for pointing that out! In that case this gives an over determined system if i x j > 2j + i. The iterative solutions would probably not converge. – russcarr Feb 20 '15 at 19:11
• But that condition implies that is over determined for pretty much all the values of $(i,j)$. If $i>3$ and $j>3$ is over determined. Check this matrix ($i$ is represented by rows and $j$ in columns) \begin{pmatrix}2 & 3 & 4 & 5 & 6 & 7 & 8 & 9\cr 2 & 2 & 2 & 2 & 2 & 2 & 2 & 2\cr 2 & 1 & 0 & -1 & -2 & -3 & -4 & -5\cr 2 & 0 & -2 & -4 & -6 & -8 & -10 & -12\cr 2 & -1 & -4 & -7 & -10 & -13 & -16 & -19\cr 2 & -2 & -6 & -10 & -14 & -18 & -22 & -26\cr 2 & -3 & -8 & -13 & -18 & -23 & -28 & -33\cr 2 & -4 & -10 & -16 & -22 & -28 & -34 & -40\end{pmatrix} – nicoguaro Feb 21 '15 at 1:26
• @nicoguaro Do you disagree? What does this matrix imply to the original system of equations? – russcarr Feb 21 '15 at 3:23
• I don't disagree. I just said that for almost every value of $i$ and $j$ that condition is true, i.e., in practical terms the system is always over constrained. The matrix is just a depiction of the inequality. – nicoguaro Feb 21 '15 at 3:32

I don't get why you need an iterative solution. Reorder to $RA =C-B$ and solve columnwise: Choose some $a_i$ and $i$th columns of $C$ and $B$ that fit. By the way: You may want to exclude the solution $A=B=C=0$.

• How do I choose those values? – nxkryptor Feb 23 '15 at 14:38
• Oops, I missed the extra structure in B and C... – Dirk Feb 27 '15 at 7:26