login  home  contents  what's new  discussion  bug reports     help  links  subscribe  changes  refresh  edit

Edit detail for numerical linear algebra revision 4 of 5

1 2 3 4 5
Editor: test1
Time: 2014/10/07 13:56:49 GMT+0
Note:

changed:
-solve(rhs(eigen.1),15)
solve(rhs(eigen.1),10.0^(-15))

changed:
-Why would you expect to be able to recover the characteristic polynomial?  There is always round off error for finite precision arithmetic. For integer approximations, it would be worse.
-The command {\tt solve: (Polynomial Fraction Integer, PositiveInteger)->List Equation Polynomial Integer}
-solves the equation over the integers, so it is {\it not} accurate. For example:
To recover characteristic polynomial one needs good precision, which means very
small second argument to solve.  Giving large second arguments means very poor
accuracy:

removed:
-ev:= solve(rhs(eigen.1),1.0*10^(-50))
-cp:= reduce(*, [rhs(x)-lhs(x) for x in ev])

added:

To have use of precise solution we also need to increase precision of other
floating point computations using digits.  Unfortunately, this leads to
ugly display, so we only show final result:
\begin{axiom}
digits(120)
ev:= solve(rhs(eigen.1),1.0*10^(-100));
cp:= reduce(*, [rhs(x)-lhs(x) for x in ev]);
rhs(eigen.1) - cp
\end{axiom}

I'm new to Axiom, so maybe I'm doing things in a stupid way.

I want to get (estimates of) the eigenvalues of a 10x10 matrix of floats:

fricas
m := matrix([[random()$Integer for i in 1..10] for j in 1..10]); sm := m + transpose(m); smf:Matrix Float := sm

\label{eq1}\left[ 
\begin{array}{cccccccccc}
{19101368.0}&{61278041.0}&{70995855.0}&{34982229.0}&{73454712.0}&{3
3749409.0}&{43151276.0}&{56873386.0}&{36952224.0}&{71646262.0}
\
{61278041.0}&{2968650.0}&{85736917.0}&{24332967.0}&{4553241.0}&{6
3473373.0}&{20558912.0}&{92297086.0}&{54375004.0}&{83468790.0}
\
{70995855.0}&{85736917.0}&{15829350.0}&{1 \<u> 19704936.0}&{365
22809.0}&{39703630.0}&{6801050.0}&{36762682.0}&{42794759.0}&{3
5393317.0}
\
{34982229.0}&{24332967.0}&{1 \</u> 19704936.0}&{26538740.0}&{489
82062.0}&{89215221.0}&{64457881.0}&{15519476.0}&{57773751.0}&{2
9448891.0}
\
{73454712.0}&{4553241.0}&{36522809.0}&{48982062.0}&{46207882.0}&{4
4673466.0}&{1 \<u> 00571933.0}&{60490393.0}&{96108269.0}&{59277
125.0}
\
{33749409.0}&{63473373.0}&{39703630.0}&{89215221.0}&{44673466.0}&{1
826916.0}&{1 \</u> 19858422.0}&{67603577.0}&{82798875.0}&{762919
26.0}
\
{43151276.0}&{20558912.0}&{6801050.0}&{64457881.0}&{1 \<u> 0057
1933.0}&{1 \</u> 19858422.0}&{1 \<u> 28173768.0}&{91100645.0}&{744
38999.0}&{68433103.0}
\
{56873386.0}&{92297086.0}&{36762682.0}&{15519476.0}&{60490393.0}&{6
7603577.0}&{91100645.0}&{53378192.0}&{36023315.0}&{94404719.0}
\
{36952224.0}&{54375004.0}&{42794759.0}&{57773751.0}&{96108269.0}&{8
2798875.0}&{74438999.0}&{36023315.0}&{1 \</u> 28828432.0}&{30785
633.0}
\
{71646262.0}&{83468790.0}&{35393317.0}&{29448891.0}&{59277125.0}&{7
6291926.0}&{68433103.0}&{94404719.0}&{30785633.0}&{91534560.0}
(1)
Type: Matrix(Float)

The problem is: If I now call eigenvalues(smf) on the symmetric float matrix smf Axiom 3.0 Beta (February 2005) runs for a very long time (uncomment code if you want to try it):

fricas
)set messages time on

Try this::

fricas
eigen:=eigenvalues(sm)

\label{eq2}\left[ \left(\%B \mid{{{\%B}^{10}}-{{514387858}\ {{\%B}^{9}}}-{{81063769717750363}\ {{\%B}^{8}}}+{{186174851844360357608590
98}\ {{\%B}^{7}}}+{{1995237067328519443299986185143419}\ {{\%B}^{6}}}-{{196489116329546029723739477630207952322692}\ {{\%B}^{5}}}-{{17713684304188911913775711272171401718653930456014}\ {{\%B}^{4}}}+{{56122609706355950813933443516197394131149130564764128
1838}\ {{\%B}^{3}}}+{{399896551509446465317296211063275390420
68724401313082575314065916}\ {{\%B}^{2}}}-{{45738107297399937
0254434039690959701937607418001547553700564621454196664}\  \%B}-{123950107976044776612898516638225104195271993434900035759953
61443410961398217808}}\right) \right](2)
Type: List(Union(Fraction(Polynomial(Integer)),SuchThat?(Symbol,Polynomial(Integer))))
fricas
Time: 0.01 (IN) + 0.01 (EV) + 0.07 (OT) = 0.09 sec
solve(rhs(eigen.1),10.0^(-15))

\label{eq3}\begin{array}{@{}l}
\displaystyle
\left[{\%B = -{14233142.7839271088 \_ 23}}, \: \right.
\
\
\displaystyle
\left.{\%B = -{54892641.0998019923 \_ 75}}, \: \right.
\
\
\displaystyle
\left.{\%B = -{93438314.6687504400 \_ 27}}, \: \right.
\
\
\displaystyle
\left.{\%B = -{104131052.5421902803 \_ 1}}, \: \right.
\
\
\displaystyle
\left.{\%B = -{139552430.7909528369 \_ 9}}, \: \right.
\
\
\displaystyle
\left.{\%B ={590433871.2921979769 \_ 2}}, \: \right.
\
\
\displaystyle
\left.{\%B ={136341415.0492967303 \_ 4}}, \: \right.
\
\
\displaystyle
\left.{\%B ={118357957.8969694609 \_ 4}}, \: \right.
\
\
\displaystyle
\left.{\%B ={51852746.7256929047 \_ 26}}, \: \right.
\
\
\displaystyle
\left.{\%B ={23649448.9214655855 \_ 95}}\right] 
(3)
Type: List(Equation(Polynomial(Float)))
fricas
Time: 0.03 (EV) + 0.01 (OT) = 0.04 sec

Thank you! This helps, but doesn't answer everything. Since interestingly:

fricas
charpol := reduce(*, [ rhs(x) - lhs(x) for x in % ])

\label{eq4}\begin{array}{@{}l}
\displaystyle
{{\%B}^{10}}-{{514387858.0}\ {{\%B}^{9}}}-{{8106376 \<u> 971775
0363.001}\ {{\%B}^{8}}}+ 
\
\
\displaystyle
{{0.1861748518 \</u> 4436035761 E 26}\ {{\%B}^{7}}}+ 
\
\
\displaystyle
{{0.1995237067 \<u> 3285194434 E 34}\ {{\%B}^{6}}}- 
\
\
\displaystyle
{{0.1964891163 \</u> 2954602973 E 42}\ {{\%B}^{5}}}- 
\
\
\displaystyle
{{0.1771368430 \<u> 4188911913 E 50}\ {{\%B}^{4}}}+ 
\
\
\displaystyle
{{0.5612260970 \</u> 6355950813 E 57}\ {{\%B}^{3}}}+ 
\
\
\displaystyle
{{0.3998965515 \<u> 0944646532 E 65}\ {{\%B}^{2}}}- 
\
\
\displaystyle
{{0.4573810729 \</u> 7399937026 E 72}\  \%B}- 
\
\
\displaystyle
{0.1239501079 \<u> 7604477661 E 80}
(4)
Type: Polynomial(Float)
fricas
Time: 0 sec

we cannot recover the characteristic polynomial from this solution: Even if a large number is passed to solve, accuracy does not increase.

To recover characteristic polynomial one needs good precision, which means very small second argument to solve. Giving large second arguments means very poor accuracy:

fricas
solve(x+11/10,3)

\label{eq5}\left[{x = - 1}\right](5)
Type: List(Equation(Polynomial(Fraction(Integer))))
fricas
Time: 0.01 (IN) = 0.01 sec

To have use of precise solution we also need to increase precision of other floating point computations using digits. Unfortunately, this leads to ugly display, so we only show final result:

fricas
digits(120)

\label{eq6}20(6)
Type: PositiveInteger?
fricas
Time: 0 sec
ev:= solve(rhs(eigen.1),1.0*10^(-100));
Type: List(Equation(Polynomial(Float)))
fricas
Time: 0.09 (EV) = 0.09 sec
cp:= reduce(*, [rhs(x)-lhs(x) for x in ev]);
Type: Polynomial(Float)
fricas
Time: 0 sec
rhs(eigen.1) - cp

\label{eq7}\begin{array}{@{}l}
\displaystyle
-{{0.5714936956 \<u> 4 E - 100}\ {{\%B}^{9}}}+ 
\
\
\displaystyle
{{0.2473701638 \</u> 8 E - 91}\ {{\%B}^{8}}}+{{0.1465027643 \<u> 5 E - 83}\ {{\%B}^{7}}}- 
\
\
\displaystyle
{{0.7490726307 \</u> 8 E - 75}\ {{\%B}^{6}}}-{{0.3974731329 E - 67}\ {{\%B}^{5}}}+ 
\
\
\displaystyle
{{0.4873818394 \<u> 7 E - 59}\ {{\%B}^{4}}}+{{0.3998713024 \</u> 3
6 E - 51}\ {{\%B}^{3}}}+ 
\
\
\displaystyle
{{0.145275416 E - 44}\ {{\%B}^{2}}}-{{0.1071484665 \<u> 3 E - 3
5}\  \%B}- 
\
\
\displaystyle
{0.1609012932 \</u> 1 E - 28}
(7)
Type: Polynomial(Float)
fricas
Time: 0 sec

fricas
A:=[[a,b],[c,d]]

\label{eq8}\left[{\left[ a , \: b \right]}, \:{\left[ c , \: d \right]}\right](8)
Type: List(List(Symbol))
fricas
Time: 0.01 (OT) = 0.01 sec

fricas
A:=matrix[[cos(x),-sin(x)],[sin(x),cos(x)]]

\label{eq9}\left[ 
\begin{array}{cc}
{\cos \left({x}\right)}& -{\sin \left({x}\right)}
\
{\sin \left({x}\right)}&{\cos \left({x}\right)}
(9)
Type: Matrix(Expression(Integer))
fricas
Time: 0.04 (OT) = 0.04 sec
A(1,1)

\label{eq10}\cos \left({x}\right)(10)
Type: Expression(Integer)
fricas
Time: 0 sec
eigen:=eigenvalues(A)
There are 1 exposed and 1 unexposed library operations named eigenvalues having 1 argument(s) but none was determined to be applicable. Use HyperDoc Browse, or issue )display op eigenvalues to learn more about the available operations. Perhaps package-calling the operation or using coercions on the arguments will allow you to apply the operation.
Cannot find a definition or applicable library operation named eigenvalues with argument type(s) Matrix(Expression(Integer))
Perhaps you should use "@" to indicate the required return type, or "$" to specify which version of the function you need.

Unfortunately, currently eigenvalues does not work for general expressions, which causes the failure above.

fricas
A:=matrix[[cos(x)-L,-sin(x)],[sin(x),cos(x)-L]]

\label{eq11}\left[ 
\begin{array}{cc}
{{\cos \left({x}\right)}- L}& -{\sin \left({x}\right)}
\
{\sin \left({x}\right)}&{{\cos \left({x}\right)}- L}
(11)
Type: Matrix(Expression(Integer))
fricas
Time: 0 sec
A(1,1)*A(2,2)

\label{eq12}{{\cos \left({x}\right)}^{2}}-{2 \  L \ {\cos \left({x}\right)}}+{{L}^{2}}(12)
Type: Expression(Integer)
fricas
Time: 0 sec
A(2,1)*A(1,2)

\label{eq13}-{{\sin \left({x}\right)}^{2}}(13)
Type: Expression(Integer)
fricas
Time: 0 sec
A(1,1)*A(2,2)-A(2,1)*A(1,2)

\label{eq14}{{\sin \left({x}\right)}^{2}}+{{\cos \left({x}\right)}^{2}}-{2 \  L \ {\cos \left({x}\right)}}+{{L}^{2}}(14)
Type: Expression(Integer)
fricas
Time: 0 sec
B := solve(A(1,1)*A(2,2)-A(2,1)*A(1,2)=0,L)

\label{eq15}\left[{L ={{{\sqrt{- 1}}\ {\sin \left({x}\right)}}+{\cos \left({x}\right)}}}, \:{L ={-{{\sqrt{- 1}}\ {\sin \left({x}\right)}}+{\cos \left({x}\right)}}}\right](15)
Type: List(Equation(Expression(Integer)))
fricas
Time: 0.01 (IN) + 0.01 (EV) + 0.01 (OT) = 0.03 sec
B(1)

\label{eq16}L ={{{\sqrt{- 1}}\ {\sin \left({x}\right)}}+{\cos \left({x}\right)}}(16)
Type: Equation(Expression(Integer))
fricas
Time: 0 sec

fricas
solve(x^2 - 2,x)

\label{eq17}\left[{{{{x}^{2}}- 2}= 0}\right](17)
Type: List(Equation(Fraction(Polynomial(Integer))))
fricas
Time: 0.01 (EV) = 0.01 sec
sqrt(2)

\label{eq18}\sqrt{2}(18)
Type: AlgebraicNumber?
fricas
Time: 0 sec
solve(x^2=4,x)

\label{eq19}\left[{x = 2}, \:{x = - 2}\right](19)
Type: List(Equation(Fraction(Polynomial(Integer))))
fricas
Time: 0.01 (IN) = 0.01 sec

fricas
P:=matrix[[a, b], [1.0 - a, 1.0 - b]]

\label{eq20}\left[ 
\begin{array}{cc}
a & b 
\
{-{{1.0}\  a}+{1.0}}&{-{{1.0}\  b}+{1.0}}
(20)
Type: Matrix(Polynomial(Float))
fricas
Time: 0 sec
eigenvectors(P)
Internal Error The function eigenvectors with signature hashcode is missing from domain EigenPackage(Float)