1
$\begingroup$

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.

UPDATE 1

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 2

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 3

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?

$\endgroup$
8
  • 1
    $\begingroup$ Why not just take $A = B = C = 0$? If you use $A = B = 0$ as your initial guesses, your iterative method will converge to a solution of $RA + B = C$ in one iteration. $\endgroup$
    – littleO
    Feb 19 '15 at 9:21
  • 1
    $\begingroup$ Are you sure this is what you have to do? Once you have chosen $A$ and $B$ you get directly a $C$ that satisfy the solution by computing $RA+B$! I don't understand why you would need an iteration. $\endgroup$ Feb 19 '15 at 9:35
  • 1
    $\begingroup$ But if the inputs are $A$, $B$ and $C$ and the unknown is $R$ then the problem becomes interesting. $\endgroup$ Feb 19 '15 at 9:37
  • 2
    $\begingroup$ @GillesBonnet For a particular row in matrix C the elements are constant irrespective of the column. $\endgroup$
    – nxkryptor
    Feb 19 '15 at 9:40
  • 1
    $\begingroup$ Is there a relation between $i$ and $j$? The structure of $A$ indicates that they are equal. If $i=j$, then you have $i^2$ equations for $3i$ unknowns. Is there a reason why this should have a nontrivial solution? $\endgroup$ Feb 19 '15 at 11:27
0
$\begingroup$

So the problem is you want to satisfy $$ R \operatorname{diag}(a) = c\mathbf{1}^T - \mathbf{1}b^T $$ for a given matrix $R$, and unknown column vectors $a,b,c$, and $\mathbf{1}$ is the vector of all $1$'s. This says that $R$ is of rank $2$.

If you know for a fact that $R$ has this form already, then you can compute its singular value decomposition $R = U\Sigma V^T$, or in dyadic form, $$ R = \sum_{i} \sigma_i u_iv_i^T $$ The singular values will be zero beyond the first two, so $$ R = \sigma_1 u_1 v_1^T + \sigma_2 u_2 v_2^T = [\operatorname{diag}(a)^{-1} c]\mathbf{1}^T - [\operatorname{diag}(a) \mathbf{1} ]b^T $$ This allows you to identify all the unknowns directly; one of your right singular vectors should all have the same elements, and the other right singular vector will be $b$. Then $a$ is the left singular vector corresponding to $b$ times minus the singular value. With the left singular vector corresponding to the constant right singular vector, and the obtained $a$, you can back out $c$.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.