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

fricas
RR ==> Expression Integer
Type: Void
fricas
CC ==> Expression Complex Integer
Type: Void
fricas
PI ==> PositiveInteger
Type: Void
fricas
NN ==> NonNegativeInteger
Type: Void
fricas
CF ==> Complex Float
Type: Void
fricas
IF ==> Interval Float
Type: Void
fricas
CIF ==> Complex Interval Float
Type: Void
fricas
VIF ==> Vector Interval Float
Type: Void
fricas
VCIF ==> Vector Complex Interval Float
Type: Void
fricas
MIF ==> Matrix Interval Float
Type: Void
fricas
free eps
Type: Void
fricas
eps:Float:=0.000000001

\label{eq1}0.1 E - 8(1)
Type: Float
fricas
free maxIter
Type: Void
fricas
maxIter:PI:=1000

\label{eq2}1000(2)
Type: PositiveInteger?
fricas
mean(x:IF):Float ==
  x0:Float:=inf(x)$IF
  x0+width(x)/2
Function declaration mean : Interval(Float) -> Float has been added to workspace.
Type: Void
fricas
iabs(x:IF):Float == abs mean(x)
Function declaration iabs : Interval(Float) -> Float has been added to workspace.
Type: Void
fricas
maxInd(k:PI,n:NN,S:MIF):PI ==
  m:PI:=k+1
  for i in k+2..n repeat
    if iabs(S(k,i)) > iabs(S(k,m)) then m:=i::PI
  return m
Function declaration maxInd : (PositiveInteger, NonNegativeInteger, Matrix(Interval(Float))) -> PositiveInteger has been added to workspace.
Type: Void
fricas
jacobi(M:MIF):Record(ev:VIF,EV:MIF) ==
  not square? M => error
  S:MIF:=copy M
  --eps:Float:=0.000000001
  n:NN:=nrows S
  i:PI; k:PI; l:PI; m:PI
  state:NN
  s:IF;c:IF;t:IF;p:IF;y:IF;d:IF;r:IF
  ind:List(PI):=[1 for i in 1..n]
  changed:List Boolean:=[false for i in 1..n]
  e:VIF:=new(n,0$IF)$VIF
  E:MIF:=new(n,n,0$IF)$MIF
  E:=diagonalMatrix [1$IF for i in 1..n]
  A:IF; B:IF
  count:NN:=0
  state:=n
  for k in 1..n repeat
    ind.k:=maxInd(k,n,S)
    e.k:=S(k,k)
    changed.k:=true
  while (state ~= 0) and (count < maxIter) repeat
    m:=1
    for k in 2..n-1 repeat
      if iabs(S(k,ind.k)) > iabs(S(m,ind.m)) then
        m:=k
    k:=m
    l:=ind.m
    p:=S(k,l)
    --
    y:=(e.l - e.k)/2
    d:=interval(iabs(y))$IF+sqrt(p^2+y^2)
    r:=sqrt(p^2+d^2)
    c:=d*recip(r)
    s:=p*recip(r)
    t:=p*p*recip(d)
    if mean(y)<0 then
      s:=-s
      t:=-t
    S(k,l):=0$IF
    --
    y:IF:=e.k
    e.k:=y-t
    if changed.k and iabs(t)<= eps then
      changed.k:=false
      state:=state-1
    else
      if (not changed.k) and iabs(t)> eps then
        changed.k:=true
        state:=state+1
    --
    y:IF:=e.l
    e.l:=y+t
    if changed.l and iabs(t)<= eps then
      changed.l:=false
      state:=state-1
    else
      if (not changed.l) and iabs(t)> eps then
        changed.l:=true
        state:=state+1
    --
    for i in 1..k-1 repeat
      A:=S(i,k); B:=S(i,l)
      S(i,k):=c*A-s*B
      S(i,l):=s*A+c*B
    for i in k+1..l-1 repeat
      A:=S(k,i); B:=S(i,l)
      S(k,i):=c*A-s*B
      S(i,l):=s*A+c*B
    for i in l+1..n repeat
      A:=S(k,i); B:=S(l,i)
      S(k,i):=c*A-s*B
      S(l,i):=s*A+c*B    
    --
    for i in 1..n repeat
      A:=E(i,k); B:=E(i,l)
      E(i,k):=c*A-s*B
      E(i,l):=s*A+c*B
    --
    ind.k := maxInd(k,n,S)
    ind.l := maxInd(l,n,S)
    --
    count:=count+1
    output([count,state])
  --
  return [e,E]$Record(ev:VIF,EV:MIF)
Function declaration jacobi : Matrix(Interval(Float)) -> Record(ev: Vector(Interval(Float)),EV: Matrix(Interval(Float))) has been added to workspace.
Type: Void
fricas
S0:= matrix [[4,-30,60,-35], [-30,300,-675,420], 
             [60,-675,1620,-1050],[-35,420,-1050,700]]

\label{eq3}\left[ 
\begin{array}{cccc}
4 & -{30}&{60}& -{35}
\
-{30}&{300}& -{675}&{420}
\
{60}& -{675}&{1620}& -{1050}
\
-{35}&{420}& -{1050}&{700}
(3)
Type: Matrix(Integer)
fricas
M:=map(s+->interval(s)$IF,S0)

\label{eq4}\left[ 
\begin{array}{cccc}
{\left[{4.0}, \:{4.0}\right]}&{\left[ -{30.0}, \: -{30.0}\right]}&{\left[{6
0.0}, \:{60.0}\right]}&{\left[ -{35.0}, \: -{35.0}\right]}
\
{\left[ -{30.0}, \: -{30.0}\right]}&{\left[{300.0}, \:{300.0}\right]}&{\left[ -{675.0}, \: -{675.0}\right]}&{\left[{420.0}, \:{420.0}\right]}
\
{\left[{60.0}, \:{60.0}\right]}&{\left[ -{675.0}, \: -{675.0}\right]}&{\left[{1
620.0}, \:{1620.0}\right]}&{\left[ -{1050.0}, \: -{1050.0}\right]}
\
{\left[ -{35.0}, \: -{35.0}\right]}&{\left[{420.0}, \:{420.0}\right]}&{\left[ -{1050.0}, \: -{1050.0}\right]}&{\left[{700.0}, \:{700.0}\right]}
(4)
Type: Matrix(Interval(Float))
fricas
R:=jacobi(M)
fricas
Compiling function mean with type Interval(Float) -> Float
fricas
Compiling function iabs with type Interval(Float) -> Float
fricas
Compiling function maxInd with type (PositiveInteger, 
      NonNegativeInteger, Matrix(Interval(Float))) -> PositiveInteger
fricas
Compiling function jacobi with type Matrix(Interval(Float)) -> 
      Record(ev: Vector(Interval(Float)),EV: Matrix(Interval(Float)))
fricas
Compiling function G710 with type Integer -> Boolean
fricas
Compiling function G712 with type NonNegativeInteger -> Boolean 
   [1, 4]
   [2, 4]
   [3, 4]
   [4, 4]
   [5, 4]
   [6, 4]
   [7, 4]
   [8, 4]
   [9, 4]
   [10, 4]
   [11, 4]
   [12, 4]
   [13, 2]
   [14, 1]
   [15, 0]

\label{eq5}\begin{array}{@{}l}
\displaystyle
\left[{
\begin{array}{@{}l}
\displaystyle
ev ={
\begin{array}{@{}l}
\displaystyle
\left[{\left[{0.1666428611 \<u> 7189038239}, \:{0.1666428611 \</u> 7189054445}\right]}, \right.
\
\
\displaystyle
\left.\:{\left[{37.1014913651 \<u> 27658067}, \:{37.1014913651 \</u> 27658308}\right]}, \: \right.
\
\
\displaystyle
\left.{\left[{2585.2538109289 \<u> 223143}, \:{2585.2538109289 \</u> 223146}\right]}, \: \right.
\
\
\displaystyle
\left.{\left[{1.4780548447 \<u> 781367581}, \:{1.4780548447 \</u> 781370694}\right]}\right] 
(5)
Type: Record(ev: Vector(Interval(Float)),EV: Matrix(Interval(Float)))
fricas
checkResult(R,M) ==
  n:=nrows M
  for i in 1..n repeat
    v:=column(R.EV,i)
    d:=M*v-R.ev.i*v
    output(sqrt(iabs(dot(d,d))))
Type: Void
fricas
-- Eigenvalues
R.ev

\label{eq6}\begin{array}{@{}l}
\displaystyle
\left[{\left[{0.1666428611 \<u> 7189038239}, \:{0.1666428611 \</u> 7189054445}\right]}, \: \right.
\
\
\displaystyle
\left.{\left[{37.1014913651 \<u> 27658067}, \:{37.1014913651 \</u> 27658308}\right]}, \: \right.
\
\
\displaystyle
\left.{\left[{2585.2538109289 \<u> 223143}, \:{2585.2538109289 \</u> 223146}\right]}, \: \right.
\
\
\displaystyle
\left.{\left[{1.4780548447 \<u> 781367581}, \:{1.4780548447 \</u> 781370694}\right]}\right] 
(6)
Type: Vector(Interval(Float))
fricas
-- Eigenvectors
R.EV

\label{eq7}\left[ 
\begin{array}{cccc}
{\left[{0.7926082911 \<u> 6404261165}, \:{0.7926082911 \</u> 64043
07336}\right]}&{\left[ -{0.1791862905 \<u> 3191911827}, \: -{0.1
791862905 \</u> 3191910664}\right]}&{\left[{0.0291933231 \<u> 7921
032378 \</u> 6}, \:{0.0291933231 \<u> 7921032385 \</u> 6}\right]}&{\left[ -{0.5820756994 \<u> 972225629}, \: -{0.5820756994 \</u> 9722221613}\right]<a class=?} \ {\left[{0.4519231209 \ 0048267041}, \:{0.4519231209 \ 00482 93329}\right]}&{\left[{0.7419177906 \ 0266856517}, \:{0.7419 177906 \ 0266860845}\right]}&{\left[ -{0.3287120558 \ 22912 41936}, \: -{0.3287120558 \ 2291241868}\right]}&{\left[{0.37 05021850 \ 6710162872}, \:{0.3705021850 \ 671018895}\right]} \ {\left[{0.3224163985 \ 8196501613}, \:{0.3224163985 \ 81965 21617}\right]}&{\left[ -{0.1002281368 \ 8300231775}, \: -{0.1 002281368 \ 8300230763}\right]}&{\left[{0.7914111458 \ 4119 456497}, \:{0.7914111458 \ 4119456652}\right]}&{\left[{0.509 5786345 \ 0180568066}, \:{0.5095786345 \ 0180598806}\right]} \ {\left[{0.2521611696 \ 8918677348}, \:{0.2521611696 \ 89186 95453}\right]}&{\left[ -{0.6382825282 \ 3465852565}, \: -{0.6 382825282 \ 3465848812}\right]}&{\left[ -{0.5145527499 \ 45 77199013}, \: -{0.5145527499 \ 4577198907}\right]}&{\left[{0.5 140482722 \ 221689925}, \:{0.5140482722 \ 2216930418}\right]} " title=" \label{eq7}\left[ \begin{array}{cccc} {\left[{0.7926082911 \ 6404261165}, \:{0.7926082911 \ 64043 07336}\right]}&{\left[ -{0.1791862905 \ 3191911827}, \: -{0.1 791862905 \ 3191910664}\right]}&{\left[{0.0291933231 \ 7921 032378 \ 6}, \:{0.0291933231 \ 7921032385 \ 6}\right]}&{\left[ -{0.5820756994 \ 972225629}, \: -{0.5820756994 \ 9722221613}\right]?} \ {\left[{0.4519231209 \ 0048267041}, \:{0.4519231209 \ 00482 93329}\right]}&{\left[{0.7419177906 \ 0266856517}, \:{0.7419 177906 \ 0266860845}\right]}&{\left[ -{0.3287120558 \ 22912 41936}, \: -{0.3287120558 \ 2291241868}\right]}&{\left[{0.37 05021850 \ 6710162872}, \:{0.3705021850 \ 671018895}\right]} \ {\left[{0.3224163985 \ 8196501613}, \:{0.3224163985 \ 81965 21617}\right]}&{\left[ -{0.1002281368 \ 8300231775}, \: -{0.1 002281368 \ 8300230763}\right]}&{\left[{0.7914111458 \ 4119 456497}, \:{0.7914111458 \ 4119456652}\right]}&{\left[{0.509 5786345 \ 0180568066}, \:{0.5095786345 \ 0180598806}\right]} \ {\left[{0.2521611696 \ 8918677348}, \:{0.2521611696 \ 89186 95453}\right]}&{\left[ -{0.6382825282 \ 3465852565}, \: -{0.6 382825282 \ 3465848812}\right]}&{\left[ -{0.5145527499 \ 45 77199013}, \: -{0.5145527499 \ 4577198907}\right]}&{\left[{0.5 140482722 \ 221689925}, \:{0.5140482722 \ 2216930418}\right]} " class="equation" src="images/8489447566440269298-16.0px.png" align="bottom" Style="vertical-align:text-bottom" width="793" height="76"/>(7)
Type: Matrix(Interval(Float))
fricas
checkResult(R,M)
fricas
Compiling function checkResult with type (Record(ev: Vector(Interval
      (Float)),EV: Matrix(Interval(Float))), Matrix(Interval(Float)))
       -> Void 
   0.5525395860_0173782996 E -10
   0.2051229717_9922636973 E -6
   0.2051229643_6060018097 E -6
   0.9692406540_0434008124 E -13
Type: Void




  Subject:   Be Bold !!
  ( 15 subscribers )  
Please rate this page: