3
$\begingroup$

Consider an equation $S(t)b(t) = a$, where $a, b(t) \in S^{n-1}$ are given and the vector $b(t)$ is continuous, i.e. its endpoint traces a continuous curve on the unit sphere. The task is to find continuous solution $S(t) \in SO(n)$ (numerically). Consider a time net $t_0 < ... < t_n$.

  1. In 2D case we can count the angle of $b(t)$ on the unit circle and then $S(t)$ will be the rotation matrix by this angle (inversed) corrected by the angle of $a$.
  2. In 3D case we can use the Euler theorem and correct the matrix $S(t_k)$ on each step of our iterational method by multiplying it by the rotation matrix $S'$ that sends $b(t_k)$ to $b(t_{k+1})$.
  3. In general case assume that $a=(1,0,...,0)^{T}$. I tried to construct some orthonormed basis ${b(t_0),e_2,...,e_n}$ at first step and to let $S(t_0) = [b(t_0),e_2,...,e_n]^{-1}$ using Gram-Schmidt process. Then at each step of our method I constructed the new one using Gram-Schmidt process and the previous basis as initial condition. It works, but it is time-consuming enough.

My question is how to obtain such $S(t)$ in general case more quickly?

$\endgroup$
  • $\begingroup$ Could you describe what are the spaces $S^{n-1}$ and $SO(n)$? $\endgroup$ – Paul Mar 21 '12 at 14:00
  • $\begingroup$ $SO(n)$ is the special orthogonal group of dimension $n$ (I'm assuming over the reals), which is the set of matrices $Q \in \mathbb{R}^{n \times n}$ such that $QQ^{T} = Q^{T}Q = I$, and $\det(Q) = 1$. $S^{n-1}$ is the $(n-1)$-dimensional unit sphere. $\endgroup$ – Geoff Oxberry Mar 21 '12 at 14:32
  • $\begingroup$ I'm curious about the description of method 3. When you say that you construct $S(t_{k+1})$ from $S(t_{k})$, are you applying Gram-Schmidt to $[b(t_{k+1}), e_{2}, \ldots, e_{n}]^{-1}$ again? What do you mean by "[using] the previous basis as initial condition"? $\endgroup$ – Geoff Oxberry Mar 21 '12 at 14:36
  • $\begingroup$ What's preventing you from using the infinitesimal generator of rotation in the general case? $\endgroup$ – Deathbreath Mar 21 '12 at 15:23
  • $\begingroup$ @GeoffOxberry, I mean that I apply the Gram-Schmidt process to $b(t_{k+1}), e_2(t_k),...,e_n(t_{k})$ and I receive the basis $b(t_{k+1}), e_2(t_{k+1}),...,e_n(t_{k+1})$. $\endgroup$ – Appliqué Mar 21 '12 at 15:53
2
$\begingroup$

Let $A(t_{k}) = [b(t_{k}), e_{2}, \ldots, e_{n}]$. If $b(t_{k+1}) = b(t_{k})$, nothing changes, so without loss of generality, assume that $b(t_{k+1}) - b(t_{k}) \neq 0$, and also assume that $A(t_{k})$ and $A(t_{k+1})$ are both invertible.

You can get $A(t_{k+1})^{-1}$ from $A(t_{k})^{-1}$ by noting that $A(t_{k+1})$ is a rank-one update of $A(t_{k})$, and using the Sherman-Morrison formula, with $u = b(t_{k+1}) - b(t_{k})$, and $v = e_{1}$; this rank-one update would be cheaper than naïvely inverting $A(t_{k+1})$ using an LU decomposition.

Then, since $A(t_{k+1})^{-1}$ is a rank-one update of $A(t_{k})^{-1}$, you can update your QR factorization (at least, that's how I would do Gram-Schmidt numerically) also using a rank-one updating scheme. One algorithm for accomplishing a rank-one update of a QR factorization can be found in Section 12.5.1 of the third edition of Matrix Computations by Golub and van Loan.

Both of these updates should reduce the complexity of $\mathcal{O}(n^{3})$ operations to $\mathcal{O}(n^{2})$ operations; you can even calculate $A(t_{0})^{-1}$ and its QR factorization in this fashion because $A(t_{0})$ is a rank-one update of the identity matrix.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Hm, $S(t_{k}) = [b(t_{k}),e_1(t_{k}),...,e_n(t_k)]^{T}$. Then I have to see only for update of decomposition $[b(t_{k+1}),e_1(t_{k}),...,e_n(t_k)] = QR$ and take $S(t_{k+1}) = Q^{T}$. Is it true? $\endgroup$ – Appliqué Mar 21 '12 at 17:47
  • $\begingroup$ That's another way you could do it. The first vector in your basis should just be scaled in a QR decomposition, so once you find the first QR decomposition (i.e., the identity matrix), you can repeatedly rank-one update it. $\endgroup$ – Geoff Oxberry Mar 21 '12 at 18:11
  • $\begingroup$ So if we receive $S(t_{k+1})$ using orthoghonalisation of $[b(t_{k+1}),e_2(t_k),...,e_n(t_k)]$ by one hand and updating a QR decomposition by the other we will get two different resulting matrices? Why the last will sent $b(t_{k+1})$ to $[1,0,...,0]^{T}$ and will be close to $S(t_k)$? I think that update algorithm doesn't garantee this $\endgroup$ – Appliqué Mar 21 '12 at 18:29
  • $\begingroup$ I don't know that the update algorithm will guarantee that the first column of the updated QR decomposition will be a scaled version of $b(t_{k+1})$, but the first update to the QR decomposition would be correct. $\endgroup$ – Geoff Oxberry Mar 21 '12 at 18:36
  • $\begingroup$ I tried the updating in MATLAB; it won't guarantee that the first column of the updated QR decomposition is a scaled version of $b(t_{k+1})$. $\endgroup$ – Geoff Oxberry Mar 21 '12 at 19:04
2
$\begingroup$

Rewrite your equation $$ b(t)=B(t)a $$

$B(t+\Delta t)=B^\prime B(t)$ where $$B^\prime=\left(\begin{matrix} b/\|b\| & \dot b/\|\dot b\| \end{matrix}\right)\left(\begin{matrix} \cos \theta & \sin \theta \\ -\sin \theta & \cos \theta \end{matrix}\right)\left(\begin{matrix} b/\|b\| & \dot b/\|\dot b\| \end{matrix}\right)^\ast +\\ I-\left(\begin{matrix} b/\|b\| & \dot b/\|\dot b\| \end{matrix}\right)\left(\begin{matrix} b/\|b\|& \dot b/\|\dot b\| \end{matrix}\right)^\ast$$ and $\theta$ such that $B^\prime$ takes $b(t)$ to $b(t+\Delta t)$. Now $S(t)=B(t)^\ast$.

B was constructed such that your coordinate frame is defined by $b$ and $\dot b$. Since all other directions are unaffected, we project onto the complement of $\mbox{span }\{b,\dot b\}$ (second term on the right) and then rotate within $b$ and $\dot b$ (first term on right).

You may also want to consider that your problem can be understood as

$$ S\dot b +\dot S b = 0$$ such that $$ S\in SO(n).$$

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ It's important to note that the differential condition is necessary, but not sufficient. $\endgroup$ – Geoff Oxberry Mar 21 '12 at 19:41
  • $\begingroup$ @GeoffOxberry: Yes, I've omitted the initial conditions as well, and of course any rotation on the complement of span $\{b,\dot b\}$ would still be a solution. Obviously $B^\prime$ is not unqiue. It would work without the projector as well. I just prefer having invariant spaces over null spaces :-) $\endgroup$ – Deathbreath Mar 21 '12 at 19:55

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.