Dear Forum,
Can somebody help me to solve the folloving matrix equation?
pps=:9!:11 NB. set print precision
pps 20
mp=: +/ . * NB. Matrix product
NB.=======================================================================================
NB. The following matrix equation have to be solve
NB. dX=(NbI-NbI*CT*NcI*C*NbI)*W+NbI*CT*NcI*Wx
NB.=======================================================================================
NB. BT = B transpose
BT=:|:B
NB. Nb = BT*B
Nb=:BT mp B
NB. NbI = Nb inverse
NbI=:%.Nb
NB. W = BT * l
W=: BT mp l
NB. C = [q0 q1 q2 q3 0]
C=:1 5$q0,q1,q2,q3,0
NB. CT = C transpose
CT=:|:C
NB. Wx = (1-(q0^2)-(q1^2)-(q2^2)-(q3^2))%2
Wx=:(1-(q0^2)-(q1^2)-(q2^2)-(q3^2))%2
NB. Nc = C * Nb * CT
Nc=:C mp Nb mp CT
NB. NcI = Nc inverse
NcI=:%.Nc
NB.=======================================================================================
]Nc=:C mp Nb mp CT
5736.3999999999996
]Nc=:(C mp Nb) mp CT
5736.3999999999996
]Nc=:C mp (Nb mp CT)
5736.3999999999996
NB.=======================================================================================
12j3":B
4.700 4.700 20.700 4.700 _3.300
_4.700 _20.700 _11.300 _27.300 _8.000
16.000 6.600 22.600 0.000 8.000
2.200 12.200 8.200 _7.800 4.200
_2.200 _18.200 6.200 _9.800 2.000
6.000 1.600 _2.400 10.000 8.000
9.200 _0.800 15.200 19.200 1.200
_9.200 _5.200 _6.800 _2.800 _8.000
6.000 _12.400 3.600 _10.000 _2.000
_10.300 _0.300 _24.300 _20.300 1.700
10.300 14.300 13.700 17.700 12.000
_14.000 6.600 _17.400 10.000 _2.000
_5.800 _15.800 _19.800 4.200 _3.800
5.800 29.800 _1.800 22.200 2.000
_14.000 _2.400 _6.400 _10.000 _12.000
$B
15 5
12j3":BT
4.700 _4.700 16.000 2.200 _2.200 6.000
9.200 _9.200 6.000 _10.300 10.300 _14.000 _5.800
5.800 _14.000
4.700 _20.700 6.600 12.200 _18.200 1.600
_0.800 _5.200 _12.400 _0.300 14.300 6.600 _15.800
29.800 _2.400
20.700 _11.300 22.600 8.200 6.200 _2.400
15.200 _6.800 3.600 _24.300 13.700 _17.400 _19.800
_1.800 _6.400
4.700 _27.300 0.000 _7.800 _9.800 10.000
19.200 _2.800 _10.000 _20.300 17.700 10.000 4.200
22.200 _10.000
_3.300 _8.000 8.000 4.200 2.000 8.000
1.200 _8.000 _2.000 1.700 12.000 _2.000 _3.800
2.000 _12.000
$BT
5 15
12j3":Nb
1222.600 623.600 1555.000 853.000 611.300
623.600 2549.800 805.200 1754.000 602.300
1555.000 805.200 2985.800 1034.600 602.300
853.000 1754.000 1034.600 2936.600 602.300
611.300 602.300 602.300 602.300 607.300
$Nb
5 5
12j6":NbI
0.004533 0.000333 _0.001772 _0.000311 _0.002826
0.000333 0.000749 _0.000128 _0.000382 _0.000573
_0.001772 _0.000128 0.001127 0.000039 0.000754
_0.000311 _0.000382 0.000039 0.000642 0.000016
_0.002826 _0.000573 0.000754 0.000016 0.004295
$NbI
5 5
12j3":W
_981.642
_1344.633
_1326.828
_1561.409
_610.603
$W
5 1
12j3":C
0.500 0.500 0.500 0.500 0.000
$C
1 5
12j3":CT
0.500
0.500
0.500
0.500
0.000
$CT
5 1
12j3":Wx
0.500
$Wx
1
12j3":Nc
5736.400
$Nc
1 1
NcI
0.00017432536085349698
$NcI
1 1
NB.=======================================================================================
NB. dX=(NbI-NbI*CT*NcI*C*NbI)*W+NbI*CT*NcI*Wx
NB. dX=(NbI-NbI mp CT mp NcI mp C mp NbI) mp W+NbI mp CT mp NcI mp Wx
NB.=======================================================================================
NB. solution step by step
dX=:(NbI-NbI mp CT mp NcI mp C mp NbI) mp W+NbI mp CT mp NcI mp Wx
NB. dX=f+i
]a=:NbI mp CT
0.0013909444659834052
0.00028633175856724745
_0.00036699466085721166
_5.9456044547029092e_6
_0.0013142105194193007
]b=:a mp NcI
2.4247689595973174e_7
4.991488713605179e_8
_6.3976476685240161e_8
_1.0364696420582438e_9
_2.291002230352313e_7
]c=:b mp C
1.2123844797986587e_7 1.2123844797986587e_7 1.2123844797986587e_7
1.2123844797986587e_7 0
2.4957443568025895e_8 2.4957443568025895e_8 2.4957443568025895e_8
2.4957443568025895e_8 0
_3.1988238342620081e_8 _3.1988238342620081e_8 _3.1988238342620081e_8
_3.1988238342620081e_8 0
_5.1823482102912192e_10 _5.1823482102912192e_10 _5.1823482102912192e_10
_5.1823482102912192e_10 0
_1.1455011151761565e_7 _1.1455011151761565e_7 _1.1455011151761565e_7
_1.1455011151761565e_7 0
]d=:c mp NbI
3.3727189656402491e_10 6.942883603207513e_11 _8.8987726198450762e_11
_1.441671712780783e_12 _3.1866568738641721e_10
6.9428836032077935e_11 1.4292217412350897e_11 _1.8318497076221252e_11
_2.9677417531211663e_13 _6.5598669749826067e_11
_8.8987726198451706e_11 _1.8318497076220693e_11 2.3479025363938921e_11
3.8037882477597397e_13 8.4078558655125846e_11
_1.4416717127807238e_12 _2.9677417531209235e_13 3.8037882477595918e_13
6.1624385209859802e_15 1.3621393066516946e_12
_3.1866568738642073e_10 _6.5598669749824141e_11 8.4078558655125924e_11
1.3621393066517597e_12 3.0108592311420753e_10
]e=:NbI - d
0.0045328344164105271 0.00033293998348083668 _0.0017724766597981138
_0.00031140912439777367 _0.0028261600905336174
0.00033293998348086812 0.00074882329681384913 _0.00012756882532661765
_0.00038153100293938686 _0.00057288268896462887
_0.0017724766597981279 _0.00012756882532660762 0.0011266169145724698
3.9439332284661476e_5 0.00075421374633090137
_0.00031140912439777405 _0.00038153100293938501 3.9439332284660466e_5
0.00064160958749499744 1.6408293152415618e_5
_0.0028261600905336603 _0.00057288268896460545 0.00075421374633090712
1.6408293152416783e_5 0.0042952959863754138
]f=:e mp W
_0.33363012653942103
_0.21892822215417146
_0.10546286817415873
_0.24545184018518121
_0.1044605740131157
]g=:NbI mp CT
0.0013909444659834052
0.00028633175856724745
_0.00036699466085721166
_5.9456044547029092e_6
_0.0013142105194193007
]h=:g mp NcI
2.4247689595973174e_7
4.991488713605179e_8
_6.3976476685240161e_8
_1.0364696420582438e_9
_2.291002230352313e_7
]i=:g mp Wx
0.00069547223299170258 0.00014316587928362372 _0.00018349733042860583
_2.9728022273514546e_6 _0.00065710525970965033
]dX=:f + i NB. step by step solution
_0.33293465430642932
_0.21878505627488784
_0.10564636550458734
_0.24545481298740857
_0.10511767927282535
NB.=======================================================================================
NB. VERSION 1
NB. dX=(NbI-NbI*CT*NcI*C*NbI)*W+NbI*CT*NcI*Wx
NB. dX=((NbI-(((NbI*CT)*NcI)*C*NbI))*W)+((NbI*CT)*NcI*Wx)
]dX=:((NbI-(((NbI mp CT) mp NcI) mp C mp NbI)) mp W)+((NbI mp CT) mp NcI mp
Wx)
_0.33363000530097303
_0.2189281971967279
_0.10546290016239707
_0.24545184070341602
_0.10446068856322722
NB.=======================================================================================
NB. VERSION 2
NB. dX=(((NbI-(((NbI*CT)*NcI)*C)*NbI)*W)+(NbI*CT)*NcI)*Wx
]dX=:(((NbI-(((NbI mp CT) mp NcI) mp C) mp NbI) mp W)+(NbI mp CT) mp NcI) mp
Wx
_0.16681494203126254 _0.10946408611964216 _0.052731466075317708
_0.12272592061082542 _0.052230401556669366
NB.=======================================================================================
NB. VERSION 3
NB. dX=((NbI-(((NbI*CT)*NcI)*C)*NbI)*W)+(NbI*CT)*NcI*Wx
]dX=:((NbI-(((NbI mp CT) mp NcI) mp C) mp NbI) mp W)+(NbI mp CT) mp NcI mp Wx
_0.33363000530097303
_0.2189281971967279
_0.10546290016239707
_0.24545184070341602
_0.10446068856322722
NB.=======================================================================================
NB. VERSION 4
NB. dX=(NbI-NbI*CT*NcI*C*NbI)*W+NbI*CT*NcI*Wx
NB. D=(NbI*CT)*NcI
NB. dX=(NbI-D*C*NbI)*W+D*Wx
NB. dX=((NbI-((D*C)*NbI))*W)+D*Wx
]dX=:(NbI-((D mp C) mp NbI)) mp W + D mp Wx[D=:(NbI mp CT) mp NcI
_0.33363012560096128
_0.21892822202521539
_0.10546286851468895
_0.24545184023593167
_0.10446057488621374
NB.=======================================================================================
NB. VERSION 5
NB. dX=(NbI-NbI*CT*NcI*C*NbI)*W+NbI*CT*NcI*Wx
NB. dX=:(NbI-NbI mp CT mp NcI mp C mp NbI) mp W+NbI mp CT mp NcI mp Wx
]dX=:(NbI-NbI mp CT mp NcI mp C mp NbI) mp W+NbI mp CT mp NcI mp Wx
_0.33363012560096128
_0.21892822202521539
_0.10546286851468895
_0.24545184023593167
_0.10446057488621374
NB.=======================================================================================
NB. ORDER OF EVALUATION J Primer tell
NB. J has a right-to-left order of evaluation.
NB. J always parenthesies from right-to-left.
NB. Remember: no verb precedence and right-to-left evaluation
NB.=======================================================================================
NB. How can solve the above matrix equation?
NB.
NB. Can somebody tell me the right solution?
NB.
NB. Thanks, Erik
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm