Re: [math-fun] Decomposing Rotation into set of sub-rotations
On 5/19/16, Warren D Smith <warren.wds@gmail.com> wrote:
it wouldn't let me post that to geometric-algebra, so you might want to post it there for me.
Or you could complain about how I'm wrong some more. I have not programmed this or tried to... but it still seems to me I'm obviously right plus was all along...
I'll take option 2, thanks. Screenshots of listing (sans comments) of program implementing algorithms under discussion at https://www.dropbox.com/s/3cr364xiw628i23/given_p0.tiff ; output of column-by-column algorithm test result at https://www.dropbox.com/s/uwsmq2habm4tfgp/given_p1.tiff ; output of triangular algorithm test result at https://www.dropbox.com/s/9v3t0nicq6huf2g/given_p2.tiff . Observe that the first result is finally diagonal, whereas the second has nonzero first column. Apologies for the excessive precision, courtesy of yet another "feature". Maple source code in attached text file, should you wish to experiment.
The procedures I gave will be stable in the sense numerical errors will grow like polynomial(N), as opposed to many procedures such as gaussian elimination with partial pivoting where error can grow exponentially(N) for unfriendly matrices.
Wilkinson famously showed gauss with complete pivoting is backwards stable in the sense of a certain bound on error growth, but I think his bound was superpolynomial although subexponential.
Good points about stability, which I shall pass on to GA list in due course. Fred Lunnon
Or you could complain about how I'm wrong some more. I have not programmed this or tried to... but it still seems to me I'm obviously right plus was all along...
FWL: I'll take option 2, thanks. Screenshots of
listing (sans comments) of program implementing algorithms under discussion at
https://www.dropbox.com/s/3cr364xiw628i23/given_p0.tiff ;
output of column-by-column algorithm test result at
https://www.dropbox.com/s/uwsmq2habm4tfgp/given_p1.tiff ;
output of triangular algorithm test result at
--ok, I immediately see from that last URL that you did it wrong. A0 to A1 looks reasonable. A1 to A2 is already wrong because note row #2 was left unaltered (should have changed) but row #1 changed (should have been unaltered). Remember we are doing adjacent Givens's as per request so it would have been only rows 2 & 3 being affected during the transition from A1 to A2 zeroing first entry in row #3. Look, I really do not want to debug your code for you line by line. But you ought to perform the most trivial inspection of your output before telling me I'm a horse. -- Warren D. Smith http://RangeVoting.org <-- add your endorsement (by clicking "endorse" as 1st step)
<< I received the attachment! How?? --rwg >> Apparently small ones can on occasion be smuggled under the wire. That program has now been substantially improved --- eg. I eventually figured out how to truncate those infernal Maple 16-digit hardware fl. pt. numbers for printing, and cleaned up the rather dodgy random orthonormal matrix generator. Source/results posted on request. On 5/19/16, Warren D Smith <warren.wds@gmail.com> wrote:
... --ok, I immediately see from that last URL that you did it wrong. A0 to A1 looks reasonable. A1 to A2 is already wrong because note row #2 was left unaltered (should have changed) but row #1 changed (should have been unaltered).
(1) Steps A0 to A2 operate on rows 3 and 4, not 1 and 2 . (2) My program executes the sequence you specified as follows --- 7 4 8 2 5 9 1 3 6 ... I did previously point out this elementary malfunction, which is the responsibility of the designer rather than the programmer. Once an A_ij to be cleared is specified, Givens' R_ij alters rows & cols i,j , intentionally or otherwise! (3) I did however summon the initiative to investigate alternative scanning sequences: NW instead of SE along diagonals, and inwards rather than outwards. Unsurprisingly, they fail as well. [In fact, I still find it rather surprising that scanning columns does work --- as you point out elsewhere, both N-ward and S-ward versions.]
Remember we are doing adjacent Givens's as per request so it would have been only rows 2 & 3 being affected during the transition from A1 to A2 zeroing first entry in row #3.
(4) Please explain how a struggling mortal is expected to a clear A_ij utilising a _single_ adjacent Givens' R_{i,i+1}(t) , unless j = i+1 . Otherwise, once away from the main sub-diagonal, we are obliged to resort to general Givens' R_ij , are we not?
Look, I really do not want to debug your code for you line by line.
A prospect which would, I expect, horrify us both equally ...
But you ought to perform the most trivial inspection of your output before telling me I'm a horse.
Nay, I look forward to others inspecting it as carefully as I have. Fred Lunnon
A minor omission from the reduction algorithm concerns the final diagonal matrix, which is not necessarily the identity matrix, but may incorporate (an even number of) -1 entries along the main diagonal. It's unclear to me whether these might somehow be absorbed into rotations earlier in the sequence; otherwise they require an extra floor(n/2) rotations to eliminate. WDS raised the question how many length m = n_C_2 sequences of Givens' rotations R_ij are "valid", reducing a general orthonormal matrix to diagonal. Denoting this number f(n) , for n = 1,...,5 , I compute f(n) = 1, 1, 3, 50, 6821; for example for n = 3 , there is a choice of [[2, 1], [3, 1], [3, 2]] , [[3, 1], [2, 1], [3, 2]] , [[3, 2], [2, 1], [3, 1]] , (restricted always to the lower triangle). If complete columns are eliminated in left-to-right order, there is an independent choice of up or down row order available for each, leading to (weak) lower bound f(n) >= 2^(n-2) . Apparently for each k , the final k_C_2 moves R_ij always have i >= n+2-k ; assuming this yields an evan more feeble upper bound, maybe f(n) <= (m!)/2 ? For n = 4 , the diagonal sequence with consecutive [i, j] = ..., [2, 1], [3, 2], [4, 3], ... never occurs; from this I think it follows that no diagonal triad occurs for any n . Fred Lunnon
On 2016-05-19 10:56, Fred Lunnon wrote:
On 5/19/16, Warren D Smith <warren.wds@gmail.com> wrote:
it wouldn't let me post that to geometric-algebra, so you might want to post it there for me.
Or you could complain about how I'm wrong some more. I have not programmed this or tried to... but it still seems to me I'm obviously right plus was all along...
I'll take option 2, thanks. Screenshots of
listing (sans comments) of program implementing algorithms under discussion at
https://www.dropbox.com/s/3cr364xiw628i23/given_p0.tiff ;
output of column-by-column algorithm test result at
https://www.dropbox.com/s/uwsmq2habm4tfgp/given_p1.tiff ;
output of triangular algorithm test result at
https://www.dropbox.com/s/9v3t0nicq6huf2g/given_p2.tiff .
Observe that the first result is finally diagonal, whereas the second has nonzero first column. Apologies for the excessive precision, courtesy of yet another "feature". Maple source code in attached text file, should you wish to experiment.
I received the attachment! How?? --rwg
The procedures I gave will be stable in the sense numerical errors will grow like polynomial(N), as opposed to many procedures such as gaussian elimination with partial pivoting where error can grow exponentially(N) for unfriendly matrices.
Wilkinson famously showed gauss with complete pivoting is backwards stable in the sense of a certain bound on error growth, but I think his bound was superpolynomial although subexponential.
Good points about stability, which I shall pass on to GA list in due course.
Fred Lunnon
participants (3)
-
Fred Lunnon -
rwg -
Warren D Smith