How to find a matrix that produced a symmetric matrix?

Question

$$M=\begin{bmatrix} m_{00} & m_{01} & m_{02}\\ m_{10} & m_{11} & m_{12} \\ m_{20} & m_{21}& m_{22} \end{bmatrix} $$

$$ M^{T}M = \begin{bmatrix} a & d & e\\ d & b & f \\ e & f& c \end{bmatrix} $$

$$ M.\begin{bmatrix} x_0 \\ y_0 \\z_0 \end{bmatrix}= \begin{bmatrix} x_1 \\ y_1 \\z_1 \end{bmatrix}$$

so it will be:

$$ a = m_{20}^{2} + m_{10}^{2} + m_{00}^{2} $$ $$ b = m_{21}^{2} + m_{11}^{2} + m_{01}^{2} $$ $$ c = m_{22}^{2} + m_{12}^{2} + m_{02}^{2} $$ $$ d = m_{20}.m_{21} + m_{10}.m_{11} + m_{00}.m_{01} $$ $$ e = m_{20}.m_{22 }+ m_{10}.m_{12} + m_{00}.m_{02} $$ $$ f = m_{21}.m_{22} + m_{11}.m_{12} + m_{01}.m_{02} $$ $$ x_1 = m_{02}.z_0 + m_{01}.y_0 + m_{00}.x_0 $$ $$ y_1 = m_{12}.z_0 + m_{11}.y_0 + m_{10}.x_0 $$ $$ z_1 = m_{22}.z_0 + m_{21}.y_0 + m_{20}.x_0 $$

Question: what is a direct solution to find $\mathbf{m_{00}, m_{01}, m_{22}}$ ?

a,b,c,d,e,f as well $\mathbf{x_0,y_0,z_0,x_1,y_1,z_1}$ are known.

I tried with iteration method. here is the python code:

# select an initial m. first one has better results
m = np.array([[ 3.0 , 0.0 , 0.0 ],
              [ 0.0 , 3.0 , 0.0 ],
              [ 0.0 , 0.0 , 3.0 ]] )
# m = np.array([[ math.sqrt(a) , 0            , 0            ] ,
#               [ 0            , math.sqrt(b) , 0            ] ,
#               [ 0            , 0            , math.sqrt(c) ] ])
trace = [[]]
m1 = m.copy()
for i in range(1000):
    # Newton's method for m00,m11,m22
    # m1[0,0] = m[0,0] - (m[2,0]**2 + m[1,0]**2 + m[0,0]**2 - a) / (2*m[0,0])
    # m1[1,1] = m[1,1] - (m[2,1]**2 + m[1,1]**2 + m[0,1]**2 - b) / (2*m[1,1])
    # m1[2,2] = m[2,2] - (m[2,2]**2 + m[1,2]**2 + m[0,2]**2 - c) / (2*m[2,2])
    # direct method for m00,m11,m22
    m1[0,0] = math.sqrt( abs( a - m[2,0]**2 - m[1,0]**2 ) )
    m1[1,1] = math.sqrt( abs( b - m[2,1]**2 - m[0,1]**2 ) )
    m1[2,2] = math.sqrt( abs( c - m[1,2]**2 - m[0,2]**2 ) )

    m1[0,1] = (d - (m[2,0]*m[2,1] + m[1,0]*m[1,1])) / m[0,0]
    m1[0,2] = (x1 - (m[0,1]*y0 + m[0,0]*x0)) / z0

    m1[1,0] = ( m[2,0]*( y1 - (m[1,2]*z0 + m[1,1]*y0) ) ) / ( z1 - m[2,2]*z0 - m[2,1]*y0 )
    m1[1,2] = (y1 - (m[1,1]*y0 + m[1,0]*x0)) / z0

    m1[2,0] = (e - (m[1,0]*m[1,2] + m[0,0]*m[0,2])) / m[2,2]
    m1[2,1] = (f - (m[1,1]*m[1,2] + m[0,1]*m[0,2])) / m[2,2]

    m = 0.01 * m1 + 0.99*m
    trace.append( m.reshape(-1).tolist() )

for a supposed Matrix: $$ M = \begin{bmatrix} 1.1 & -0.1 & 0.2\\ 0.3 & 1.2 & -0.1 \\ 0.1 & -0.3 & 0.9 \end{bmatrix} $$

the result is: $$ M = \begin{bmatrix} 1.14312 & 0.209234 & 0.2\\ 0.00182697 & 1.16508 & -0.1 \\ 0.0572922 & -0.372589 & 0.9 \end{bmatrix} $$

here you can see trace of results

Answer 1

Here are two different approaches 1) and 2).

• 1. Let $A$ be a real symmetric matrix. We are looking for a square matrix $M$ of the same size as $A$ with real entries such that

$$M^TM=A$$

A necessary condition is that symmetric matrix $A$ must have $\ge 0$ eigenvalues (see explanation below).

In this case, we can diagonalize $A$ under the form

$$A=Q^T\Lambda Q\tag{1}$$

where $\Lambda$ is a diagonal matrix with $\lambda_k \ge 0$ entries. Therefore $\Lambda = \Delta^T\Delta$ where $\Delta$ is a diagonal matrix with (diagonal) entries $\sqrt{\lambda_k}$

Therefore (1) can be written

$$A=(\Delta Q)^T(\Delta Q)$$

meaning that $M=\Delta Q$ is a solution.

More generaly, for any orthogonal matrix $\Omega$,

$$M=\Omega\Delta Q$$

is a solution.

See https://math.stackexchange.com/a/2858336/305862 for a similar answer to a different question adressing rectangular matrices with rank condition.

---

Explanation about the positivity of eigenvalues :

Quadratic form $X^T(M^TM)X=(MX)^T(MX)=\|MX\|^2 \ge 0$

is (indeed!) semi-definite positive.

---

• 2. A completely different approach is to use Cholesky decomposition $A=C^TC$ (or in case one or more of the eigenvalues is $0$, the incomplete Cholesky factorization) with $C$ an upper triangular matrix.