Sylvester

Quaternionic Sylvester Equation for the HP-41

Overview

1°) Program#1
2°) Program#2 ( faster )
3°) Program#3

-These programs solve the linear equation: a q + q b = c where a , b , c are 3 given quaternions and q is unknown.

1°) Program#1

-The 1st routine solves the real 4x4 linear system M.q = c where

               a₁+b₁   -a₂-b₂   -a₃-b₃    -a₄-b₄
   M =    a₂+b₂    a₁+b₁   -a₄+b₄    a₃-b₃
               a₃+b₃     a₄-b₄    a₁+b₁   -a₂+b₂
               a₄+b₄   -a₃+b₃    a₂-b₂     a₁+b₁

with a = a₁ + a₂ i + a₃ j + a₄ k & b = a₁ + b₂ i + b₃ j + b₄ k

Data Registers: R00 = Det ( Registers R01 thru R12 are to be initialized before executing "SYLV" )

                                      • R01 = a₁              • R05 = b₁              • R09 = c₁              R13 thru R28: temp
                                      • R02 = a₂              • R06 = b₂              • R10 = c₂
                                      • R03 = a₃              • R07 = b₃              • R11 = c₃
                                      • R04 = a₄              • R08 = b₄              • R12 = c₄
Flags: /
Subroutine:    "LS3" ( cf "Linear & Non-Linear Systems for the HP-41" )

01 LBL "SYLV"
02 RCL 01
03 RCL 05
04 +
05 STO 13
06 STO 18
07 STO 23
08 STO 28
09 RCL 02
10 RCL 06
11 +
12 STO 14

13 CHS
14 STO 17
15 RCL 03
16 RCL 07
17 +
18 STO 15
19 CHS
20 STO 21
21 RCL 04
22 RCL 08
23 +
24 STO 16

25 CHS
26 STO 25
27 RCL 04
28 LASTX
29 -
30 STO 19
31 CHS
32 STO 22
33 RCL 03
34 RCL 07
35 -
36 STO 26

37 CHS
38 STO 20
39 RCL 02
40 RCL 06
41 -
42 STO 24
43 CHS
44 STO 27
45 9.00102
46 REGMOVE
47 4
48 5

49 XEQ "LS3"
50 X=0?
51 " ?"
52 X=0?
53 AVIEW
54 RCL 04
55 RCL 03
56 RCL 02
57 RCL 01
58 END

( 94 bytes / SIZE 029 )

STACK	INPUTS	OUTPUTS
T	/	t
Z	/	z
Y	/	y
X	/	x

With q = x + y i + z j + t k

Example: Solve ( 2 - 3 i + 4 j - 7 k ) q + q ( 3 + 4 i - 5 j + 6 k ) = 1 + 2 i - 3 j + 4 k

   2   STO 01     3   STO 05     1   STO 09
-3 STO 02     4   STO 06     2   STO 10
   4   STO 03   -5   STO 07    -3   STO 11
-7 STO 04     6   STO 08     4   STO 12

XEQ "SYLV"    >>>>   0.239980450                  ---Execution time = 30s---
                           RDN    0.418866080
                           RDN   -0.576246334
                           RDN    0.795210166

-Thus, q = 0.239980450 + 0.418866080 i - 0.576246334 j + 0.795210166 k

R00 = Det = 8184

Note:

-If R00 = 0, the HP-41 displays " ?"
-There is no solution or an infinite set of solutions.

2°) Program#2

-This variant is much faster and only uses quaternionic algebra. The solution is given by:

q = ( 2 Ré(b) + a + | b |² a^-1 )^-1 ( c + a^-1 c Conj(b) )

Data Registers: R00: temp ( Registers R01 thru R12 are to be initialized before executing "SYLV" )

                                      • R01 = a₁              • R05 = b₁              • R09 = c₁              R13 thru R20: temp
                                      • R02 = a₂              • R06 = b₂              • R10 = c₂
                                      • R03 = a₃              • R07 = b₃              • R11 = c₃
                                      • R04 = a₄              • R08 = b₄              • R12 = c₄
Flags: /
Subroutines:    "1/Q" & "Q*Q"   ( cf "Quaternions for the HP-41" )

01 LBL "SYLV"
02 RCL 05
03 X^2
04 RCL 06
05 X^2
06 RCL 07
07 X^2
08 RCL 08
09 X^2
10 +
11 +
12 +
13 STO 13
14 RCL 09
15 STO 17
16 RCL 10
17 STO 18
18 RCL 11
19 STO 19

20 RCL 12
21 STO 20
22 RCL 04
23 RCL 03
24 RCL 02
25 RCL 01
26 XEQ "1/Q"
27 STO 09
28 RDN
29 STO 10
30 RDN
31 STO 11
32 RDN
33 STO 12
34 RDN
35 X<> 13
36 ST* 13
37 ST* T
38 ST* Z

39 *
40 ST+ 02
41 RDN
42 ST+ 03
43 X<>Y
44 ST+ 04
45 RCL 13
46 RCL 05
47 ST+ X
48 +
49 ST+ 01
50 RCL 04
51 RCL 03
52 RCL 02
53 RCL 01
54 XEQ "1/Q"
55 STO 13
56 RDN
57 STO 14

58 RDN
59 STO 15
60 X<>Y
61 STO 16
62 1
63 CHS
64 ST* 06
65 ST* 07
66 ST* 08
67 9.001004
68 STO 00
69 8
70 +
71 REGMOVE
72 XEQ "Q*Q"
73 STO 05
74 RDN
75 STO 06
76 RDN

77 STO 07
78 X<>Y
79 STO 08
80 RCL 00
81 REGMOVE
82 XEQ "Q*Q"
83 ST+ 17
84 RDN
85 ST+ 18
86 RDN
87 ST+ 19
88 X<>Y
89 ST+ 20
90 13.001008
91 REGMOVE
92 XEQ "Q*Q"
93 END

( 161 bytes / SIZE 021 )

STACK	INPUTS	OUTPUTS
T	/	t
Z	/	z
Y	/	y
X	/	x

With q = x + y i + z j + t k

Example: Solve ( 2 - 3 i + 4 j - 7 k ) q + q ( 3 + 4 i - 5 j + 6 k ) = 1 + 2 i - 3 j + 4 k

   2   STO 01     3   STO 05     1   STO 09
-3 STO 02     4   STO 06     2   STO 10
   4   STO 03   -5   STO 07    -3   STO 11
-7 STO 04     6   STO 08     4   STO 12

XEQ "SYLV"    >>>>   0.239980450                  ---Execution time = 10s---
                           RDN    0.418866080
                           RDN   -0.576246334
                           RDN    0.795210166

-Thus, q = 0.239980450 + 0.418866080 i - 0.576246334 j + 0.795210166 k

Notes:

-Of course, this method cannot be used if a = 0
-But if | a | << | b | , it could be preferable to use another equivalent formula to avoid a loss of accuracy.
-The 3rd routine does that:

3°) Program#3

-Actually, the solution may be computed by 2 formulae:

• q = ( 2 Ré(b) + a + | b |² a^-1 )^-1 ( c + a^-1 c Conj(b) ) ( which was employed in the second program )

• q = ( c + Conj(a) c b^-1 ) ( 2 Re(a) + b + | a |² b^-1 )^-1

-To get a better accuracy, the 1st formula is used if | b | <= | a |
and the second one if | b | > | a |

Data Registers: R00: temp ( Registers R01 thru R12 are to be initialized before executing "SYLV" )

                                      • R01 = a₁              • R05 = b₁              • R09 = c₁              R13 thru R21: temp
                                      • R02 = a₂              • R06 = b₂              • R10 = c₂
                                      • R03 = a₃              • R07 = b₃              • R11 = c₃
                                      • R04 = a₄              • R08 = b₄              • R12 = c₄
Flag: F00 is cleared
Subroutines:    "1/Q" & "Q*Q"   ( cf "Quaternions for the HP-41" )

01 LBL "SYLV"
02 CF 00
03 RCL 01
04 X^2
05 RCL 02
06 X^2
07 RCL 03
08 X^2
09 RCL 04
10 X^2
11 +
12 +
13 +
14 RCL 05
15 X^2
16 RCL 06
17 X^2
18 +
19 RCL 07
20 X^2
21 +
22 RCL 08
23 X^2
24 +
25 X>Y?

26 SF 00
27 X>Y?
28 X<>Y
29 STO 13
30 1.005004
31 STO 00
32 FS? 00
33 REGSWAP
34 RCL 09
35 STO 17
36 RCL 10
37 STO 18
38 RCL 11
39 STO 19
40 RCL 12
41 STO 20
42 RCL 04
43 RCL 03
44 RCL 02
45 RCL 01
46 XEQ "1/Q"
47 STO 09
48 RDN
49 STO 10
50 RDN

51 STO 11
52 RDN
53 STO 12
54 RDN
55 X<> 13
56 ST* 13
57 ST* T
58 ST* Z
59 *
60 ST+ 02
61 RDN
62 ST+ 03
63 X<>Y
64 ST+ 04
65 RCL 13
66 RCL 05
67 ST+ X
68 +
69 ST+ 01
70 RCL 04
71 RCL 03
72 RCL 02
73 RCL 01
74 XEQ "1/Q"
75 STO 13

76 RDN
77 STO 14
78 RDN
79 STO 15
80 X<>Y
81 STO 16
82 1
83 CHS
84 ST* 06
85 ST* 07
86 ST* 08
87 9.001004
88 STO 21
89 8
90 +
91 REGMOVE
92 RCL 00
93 FS? 00
94 REGSWAP
95 XEQ "Q*Q"
96 STO 05
97 RDN
98 STO 06
99 RDN
100 STO 07

101 X<>Y
102 STO 08
103 RCL 21
104 REGMOVE
105 RCL 00
106 FS? 00
107 REGSWAP
108 XEQ "Q*Q"
109 ST+ 17
110 RDN
111 ST+ 18
112 RDN
113 ST+ 19
114 X<>Y
115 ST+ 20
116 13.001008
117 REGMOVE
118 RCL 00
119 FS?C 00
120 REGSWAP
121 XEQ "Q*Q"
122 END

( 209 bytes / SIZE 022 )

STACK	INPUTS	OUTPUTS
T	/	t
Z	/	z
Y	/	y
X	/	x

With q = x + y i + z j + t k

Example: Solve ( 2 - 3 i + 4 j - 7 k ) q + q ( 3 + 4 i - 5 j + 6 k ) = 1 + 2 i - 3 j + 4 k

   2   STO 01     3   STO 05     1   STO 09
-3 STO 02     4   STO 06     2   STO 10
   4   STO 03   -5   STO 07    -3   STO 11
-7 STO 04     6   STO 08     4   STO 12

XEQ "SYLV"    >>>>   0.239980450                  ---Execution time = 11-12s---
                           RDN    0.418866080
                           RDN   -0.576246334
                           RDN    0.795210166

-And, q = 0.239980450 + 0.418866080 i - 0.576246334 j + 0.795210166 k

References:

[1] Georgi Georgiev, Ivan Ivanov, Milena Mihaylova, Tsvetelina Dinkova - "An algorithm for solving a Sylvester quaternion equation"
[2] Drahoslava Janovsk´a, Gerhard Opfer - "Linear equations in quaternionic variables"

hp41programs

Quaternionic Sylvester Equation for the HP-41