1. Matrix Algebra
- Posted by Arthur Adamson <euclid at HOST.HORANDATA.NET>
Dec 07, 1996
-
Last edited Dec 08, 1996
--=====================_818393689==_
From Arthur Adamson re Matrix Algebra
--=====================_818393689==_
2. Re: Matrix Algebra
- Posted by Robert Craig <72614.1667 at COMPUSERVE.COM>
Dec 08, 1996
-
Last edited Dec 09, 1996
Arthur Adamson writes:
> global type vector(integer len)
> --type vector(x) creates and initializez ax x element vector to 0
>
> sequence s
> s = {0.0}
>
> for i = 2 to len by 1 do s = s & {0}
> end for
> return s
> end type
In Euphoria a "type" should return 0 (FALSE) or non-zero (TRUE).
Perhaps you mean "function".
In this code and some similar code that you posted I would rather
see something like:
repeat(0, x) -- standard Euphoria library routine, repeats 0 x times
rather than a "&" or append() loop. repeat() is simpler, executes faster, and
handles the case of length 0 correctly. (Of course by initializing s
to the empty sequence, {}, and running the for-loop from 1 to n,
the loop method above could also handle length-0 correctly).
You can initialize a "2-dimensional" sequence using:
z = repeat(repeat(0,x), y)
and so on for 3-d, 4-d etc.
Other than that, I think a matrix algebra package could be quite useful.
I don't know of anyone who is working on one.
Regards,
Rob Craig
Rapid Deployment Software
3. Re: Matrix Algebra
At 09:08 PM 12/8/96 EST, you wrote:
>---------------------- Information from the mail header -----------------------
>Sender: Euphoria Programming for MS-DOS <EUPHORIA at
>MIAMIU.ACS.MUOHIO.EDU>
>Poster: Robert Craig <72614.1667 at COMPUSERVE.COM>
>Subject: Re: Matrix Algebra
>-------------------------------------------------------------------------------
>
>Arthur Adamson writes:
>> global type vector(integer len)
>> --type vector(x) creates and initializez ax x element vector to 0
>>
>> sequence s
>> s = {0.0}
>>
>> for i = 2 to len by 1 do s = s & {0}
>> end for
>> return s
>> end type
>
>In Euphoria a "type" should return 0 (FALSE) or non-zero (TRUE).
>Perhaps you mean "function".
>
>In this code and some similar code that you posted I would rather
>see something like:
>
> repeat(0, x) -- standard Euphoria library routine, repeats 0 x times
>
>rather than a "&" or append() loop. repeat() is simpler, executes faster, and
>handles the case of length 0 correctly. (Of course by initializing s
>to the empty sequence, {}, and running the for-loop from 1 to n,
>the loop method above could also handle length-0 correctly).
>
>You can initialize a "2-dimensional" sequence using:
>
> z = repeat(repeat(0,x), y)
>
>and so on for 3-d, 4-d etc.
>
>Other than that, I think a matrix algebra package could be quite useful.
>I don't know of anyone who is working on one.
>
>Regards,
> Rob Craig
> Rapid Deployment Software
>
-------------------------------
Rob, thanks for the help. I will incorporate your suggestions but not
real soon since I will be gone for a month or so. Bye Art
Arthur P. Adamson, The Engine Man, euclid at mail.horandata.net
4. Re: Matrix Algebra
I needed to quickly put together Vector Add/Subtract/Dot/X product.
Created the following but it's not error free; I am still working on=20
ironing out some details; the arccos/arcsin would probably be better=20
done through an external call to 'C' function, I'll wait for that
feature. Having done some work in APL, I think there is alot of value
to euphoria in matrix manipulations; I was impressed by my ability to
quickly generate a 64 billion 3D array of integers, something I would be
hardpressed to do in 'C'. The special ascii characters are: 224,225,231
for
alpha, beta, gamma unit angles on the vectors.
-WWW (that's not world wide web, just my initials)
include get.e
include graphics.e
include matx.e
--Calculate magnitude of 3d vector
function SizeVect(sequence v)
atom siz
sequence r1
=20
r1 =3D power(v,2)
siz =3D power(r1[1] + r1[2] + r1[3],.5)
return siz
end function
--Print sequence vector in ijk format (with highlight)
procedure PrintVect(sequence v, integer r, integer c)
sequence sign
sign =3D " + "
position(r,c)
printf(1,"%g",v[1])
text_color(BRIGHT_CYAN)
puts(1,"i")
text_color(WHITE)
if v[2] < 0 then
sign =3D " - "
v[2] =3D v[2] * -1
end if
printf(1,sign & "%g",v[2])
text_color(BRIGHT_CYAN)
puts(1,"j")
text_color(WHITE)
sign =3D " + "
if v[3] < 0 then
sign =3D " - "
v[3] =3D v[3] * -1
end if
=20
printf(1,sign & "%g",v[3])
text_color(BRIGHT_CYAN)
puts(1,"k")
text_color(WHITE)
end procedure
--Provide data retrieval
procedure Xtern(atom in)
integer fn
fn =3D open("tst.in","w")
printf(fn,"%g",in)
close(fn)
end procedure
sequence ans, mat, mat1, mat2, wrk
atom NOTDONE, rc, scalar, w3, w4
atom t
NOTDONE =3D 1
while NOTDONE do
clear_screen()
position(20,1) =20
puts(1,"Enter 3D vector A sequence: {i,j,k} 0 to end->")
ans =3D get(0)
if atom(ans[2]) then
exit
end if
rc =3D length(ans[2])
if rc =3D 0 then
NOTDONE =3D 0
else
mat1 =3D ans[2]
w3 =3D SizeVect(mat1)
position(2,1)
puts(1,"Vector A is: ")
PrintVect(mat1,2,15)
printf(1," mag=3D%g",w3)
printf(1," =E0=3D%g =E1=3D%g =E7=3D%g",{arccos(mat1[1]/w3),arcsin(mat1[2]=
/w3),
arccos(mat1[3]/w3)})
position(21,1)
puts(1,"Enter 3D vector B sequence ->")
ans =3D get(0)
mat2 =3D ans[2]
w4 =3D SizeVect(mat2)
position(3,1)
puts(1,"Vector B is: ")
PrintVect(mat2,3,15)
printf(1," mag=3D%g",w4)
printf(1," =E0=3D%g =E1=3D%g =E7=3D%g",{arccos(mat2[1]/w4),arcsin(mat2[2]=
/w4),
arccos(mat2[3]/w4)})
-- Vector manipulations -----
-- Vector Add
mat =3D mat1 + mat2
position(5,1)
puts(1,"^A + ^B =3D ")
PrintVect(mat,5,15)
=20
-- Vector Subtract -------
mat =3D mat1 - mat2
position(6,1)
puts(1,"^A - ^B =3D ")
PrintVect(mat,6,15)
-- Dot Product ------
mat =3D mat1 * mat2
scalar =3D mat[1] + mat[2] + mat[3]
position(7,1)
puts(1,"^A =F8 ^B =3D ") -- graphic dot is 248 decimal
printf(1,"%g",scalar)
-- Cross product ----
wrk =3D repeat(0,6) -- initialize work matrix
wrk[1] =3D mat1[1] * mat2[2]
wrk[2] =3D mat1[1] * mat2[3] * -1
wrk[3] =3D mat1[2] * mat2[1] * -1
wrk[4] =3D mat1[2] * mat2[3]
wrk[5] =3D mat1[3] * mat2[1]
wrk[6] =3D mat1[3] * mat2[2] * -1
mat[1] =3D wrk[4] + wrk[6] -- new i
mat[2] =3D wrk[2] + wrk[5] -- new j
mat[3] =3D wrk[1] + wrk[3] -- new k
position(8,1)
puts(1,"^A X ^B =3D ")
PrintVect(mat,8,15)
t =3D scalar / (w3 * w4)
Xtern(t)
position(9,1)
printf(1,"Angle =E9 =3D %g degrees", arccos(scalar / (w3 * w4)))
rc =3D wait_key()
end if
end while
clear_screen()
puts(1,"--The End--")
---- Include file matx.e -----
global constant PI =3D 3.14159265
--constant PI2 =3D 1.570796326794896619231322 -- pi / 2
--constant PI3 =3D 1.047197551196597746154214 -- pi / 3
--constant PI4 =3D 0.7853981633974483096156608 -- pi / 4
--constant SQRTPI =3D 1.772453850905516027298167 -- sqrt(pi)
--constant SQRT2PI =3D 2.506628274631000502415765 -- sqrt(2.pi)
--constant TWOPI =3D 6.283185307179586476925287 -- 2.pi
--constant LNPI =3D 1.144729885849400174143427 -- ln(pi)
--constant LOGPI =3D 0.4971498726941338543512683 -- log(pi)
--constant LOGE =3D 0.4342944819032518276511289 -- log(e)
--constant LN10 =3D 2.302585092994045684017991 -- ln(10)
global constant E =3D 2.718281828459045235360287 -- exp(1)
global constant ONERAD =3D 57.295779513082320876798155 -- 1 rad in =F8
global constant ONEDEG =3D 0.017453292519943295769237 -- 1=F8 in rad
sequence cos_angls
sequence sin_angls
atom cos_init
atom sin_init
cos_init =3D 0
sin_init =3D 0
--Generate an integer sequence
-- behaves as APL interval function
global function interval(atom lim)
sequence n
n =3D {}
puts(1,"Sequence initialized\n")
for x=3D1 to lim do
n =3D append(n,x)
end for
return n
end function
global function stripCR(sequence in)
sequence wrk
wrk =3D in[1..(length(in)-1)]
return wrk
end function
--arc trig functions do a onetime build of sin/cos table
-- then extrapolates to return angle in degrees=20
-- compute arccos, return angle in degrees whose cosine is theta
global function arccos(atom in)=20
atom sel, wrk1, wrk2, wrk3
if cos_init =3D 0 then
cos_init =3D 9
cos_angls =3D repeat(0.0,360)
for x=3D1 to 360 do
cos_angls[x] =3D cos(x * ONEDEG)
if x =3D 90 or x =3D 270 then
cos_angls[x] =3D 0
end if
end for
end if
sel =3D 1=20
while in < cos_angls[sel] do
sel =3D sel + 1
end while=20
sel =3D sel - 1
wrk1 =3D in - cos_angls[sel]
wrk2 =3D cos_angls[sel+1]-cos_angls[sel]
wrk3 =3D sel + wrk1 / wrk2
return wrk3
end function
-- compute arcsin, return angle in degrees whose sine is theta
global function arcsin(atom in)=20
atom sel, wrk1, wrk2, wrk3
if sin_init =3D 0 then
sin_init =3D 9
sin_angls =3D repeat(0.0,360)
for x=3D1 to 360 do
sin_angls[x] =3D sin(x * ONEDEG)
if x =3D 90 or x =3D 360 then
sin_angls[x] =3D 0
end if
end for
end if
sel =3D 1=20
while in > sin_angls[sel] do
sel =3D sel + 1
if sin_angls[sel] < 0 then
in =3D -in
else
if sin_angls[sel-1] <=3D 0 then
in =3D -in
end if =20
end if
end while=20
sel =3D sel - 1
if sel =3D 0 then
sel =3D 359
end if
wrk1 =3D in - sin_angls[sel]
wrk2 =3D sin_angls[sel+1]-sin_angls[sel]
wrk3 =3D sel + wrk1 / wrk2
return wrk3
end function
