# Finding multiple solutions to $a = Ra$ with Symmetric Orthogonal $R \in SO(3), R = R^T$ and $a \in Re(3)$

I have a problem where I wish to find multiple solutions to the equation:

$$\mathbf{a} = \mathbf{R} \mathbf{a}$$

where I wish to calculate the solutions for $$\mathbf{R}$$ that is both symmetric ($$\mathbf{R} = \mathbf{R}^T$$) and special orthogonal ($$\mathbf{R} \in SO(3)$$) from a given $$\mathbf{a}$$ that is a real-valued 3-dimensional vector of non-zero magnitude.

For example, if I have the following vector:

$$\mathbf{a} = \begin{bmatrix} 1.70349 && -1.67761 && 9.51419\end{bmatrix}^T$$

Then two solutions are $$\mathbf{R} = \mathbf{I}$$ and

$$\mathbf{R} = \begin{bmatrix} -0.939693 && -0.0593912 && 0.336824 \\ -0.0593912 && -0.941511 && -0.331707 \\ 0.336824 && -0.331707 && 0.881204 \end{bmatrix}$$

This second solution can be generated from $$\mathbf{R} = \mathbf{Q} \mathbf{R}_z \mathbf{Q}^T$$ where $$\mathbf{R}_z$$ is a rotation about the z-axis by $$\pi$$ and

$$\mathbf{Q} = \begin{bmatrix} 0.984808 && 0.0301537 && -0.17101 \\ 0 && 0.984808 && 0.173648 \\ 0.173648 && -0.17101 && 0.969846 \end{bmatrix}$$

How would one compute the non-trivial solution? Does this approach generalize to an extended problem where the vectors are different, i.e. solving for the solution(s) for $$\mathbf{R}$$ in the equation $$\mathbf{a} = \mathbf{R} \mathbf{b}$$, given $$\mathbf{a}$$ and $$\mathbf{b}$$ ?

(The following will be contingent on $$a$$ being a unit vector, but you can easily get to that case by normalizing.)

In an orthonormal basis $$B$$ containing $$a$$, this is easy. Let's say $$a$$ is the third basis vector, and that its orientation is the same as of the standard basis (so the basis change matrix is in $$\operatorname{SO}(3)$$, though that's only cosmetic), then the matrix representation $$\tilde R$$ wrt this basis must be of the form

$$\begin{pmatrix}\cos\varphi&-\sin\varphi&0\\ \sin\varphi&\cos\varphi&0\\ 0&0&1\end{pmatrix}$$

for some $$\varphi\in[0,2\pi)$$. The condition $${\tilde R}^T=\tilde R$$ enforces $$\varphi\in\{0,\pi\}$$, one of which results in the identity matrix, while the other result in the solution you described after a change of basis: Let $$Q$$ be the change of basis matrix from the standard basis to our orthonormal basis $$B$$. Then

$$R=Q^T\tilde RQ.$$

The case $$a=Rb$$ is harder. I was thinking of considering the matrix representation wrt two different orthonormal bases, one containing $$a$$, the other $$b$$, and then doing the same as above. But this way, $$R^T=R$$ might no longer be true after a change of bases. So I'm out of ideas for now.

A nice approach is to use the fact that every $$n\times n$$ real orthogonal matrix is generated by the product of at most $$n$$ reflection (Householder) matrices.
(all vectors have unit length )
i.e. $$Q^{(k)}=I-2\mathbf v_k\mathbf v_k^T$$
$$\det\big(Q^{(k)}\big) = -1$$

In your problem $$n=3$$ and $$\det\big(R\big) = 1$$ so $$R$$ must be the product of exactly 2 Householder matrices. (Note the involution $$\big(Q^{(k)}\big)^2=I$$ so even the identity matrix may be written as a product of 2 Householder matrices.)

Select the mutually orthonormal set $$\{\mathbf a, \mathbf v_1,\mathbf v_2\}$$ so
$$\Big(R\Big)\mathbf a=\Big(\big(I-2\mathbf v_1\mathbf v_1^T\big)\big(I-2\mathbf v_2\mathbf v_2^T\big)\Big)\mathbf a=\big(I-2\mathbf v_1\mathbf v_1^T\big)\mathbf a=\mathbf a$$

You have an incremental constraint of $$R=R^T$$ (or equivalently the involution $$R^2=I$$). This is satisfied in the above since we ensured $$\mathbf v_1^T\mathbf v_2=0$$

$$R$$
$$=\big(I-2\mathbf v_1\mathbf v_1^T\big)\big(I-2\mathbf v_2\mathbf v_2^T\big)$$
$$= I -2\mathbf v_1\mathbf v_1^T - 2\mathbf v_2\mathbf v_2^T$$
$$= \big(I-2\mathbf v_2\mathbf v_2^T\big)\big(I-2\mathbf v_1\mathbf v_1^T\big)$$
$$=R^T$$

as for the second question
$$R \mathbf{b}= \mathbf{a}$$
First select the Householder matrix such that
$$\big(I-2\mathbf v_2\mathbf v_2^T\big)\mathbf b=\mathbf a$$
(check: $$\mathbf v_2 \propto \mathbf a -\mathbf b$$)

then select $$\mathbf v_1$$ such that $$\mathbf v_1 \perp \mathbf a$$ and $$\mathbf v_1 \perp \mathbf v_2$$,

Thus
$$\big(I-2\mathbf v_1\mathbf v_1^T\big)\big(I-2\mathbf v_2\mathbf v_2^T\big)\mathbf b=\big(I-2\mathbf v_1\mathbf v_1^T\big)\Big(\big(I-2\mathbf v_2\mathbf v_2^T\big)\mathbf b\Big)=\big(I-2\mathbf v_1\mathbf v_1^T\big)\mathbf a = \mathbf a$$

so
$$R:=\big(I-2\mathbf v_1\mathbf v_1^T\big)\big(I-2\mathbf v_2\mathbf v_2^T\big)$$ and
$$R=R^T$$
as desired

The answer can be found using a variation of Wahba's Problem.

Since R is symmetric, it implies that:

\begin{align} a &= Rb \\ b &= Ra \end{align}

Following Ref [1], we may then construct a symmetric matrix:

$$B = \frac{1}{2} \left(ab^T + ba^T \right)$$

And compute the SVD:

$$B = U S V^T$$

Since $$B$$ is symmetric, $$U = V$$.

If $$B$$ is of rank 2 (i.e. $$a \neq b$$), then we may compute the solution as:

$$R = U_+ V_+^T$$

where $$U_+$$ and $$V_+$$ are forced to be special orthogonal by:

\begin{align} U_+ &= U \operatorname{diag} \left( \begin{bmatrix} 1 && 1 && \det(U)\end{bmatrix} \right) \\ V_+ &= V \operatorname{diag} \left( \begin{bmatrix} 1 && 1 && \det(V)\end{bmatrix} \right) \\ \end{align}

If $$B$$ is of rank 1 (i.e. $$a=b$$) then we may introduce a single axis rotation $$W$$ about the x-axis:

$$R = U_+ W V_+^T$$

To keep R symmmetric, we have two choices of W:

\begin{align} W &= I \\ W &= \operatorname{diag} \left( \begin{bmatrix} 1 && -1 && -1 \end{bmatrix} \right) \end{align}

Which corresponds as the two variants that we seek.

References: [1] Markley, F. L. Attitude Determination using Vector Observations and the Singular Value Decomposition, Journal of the Astronautical Sciences, 1988, 38:245-258