/usr/share/euler/progs/x.e is in euler 1.61.0-10.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 | comment
Exact solver
endcomment
.. ### Exact arithmetic ###
function xlgs (A,b,n=20)
## Compute the solution of Ax=b.
## n is the number of residual iterations.
## The iteration stops, when there is no more change in x.
## You can specify an epsilon eps with eps=... as last parameter.
{LU,rows,c,d}=lu(A); LU=LU[rows];
if (isvar("eps")); localepsilon(eps); endif;
v=lusolve(LU,b[rows]);
loop 1 to n;
r=residuum(A,v,b);
vn=v-lusolve(LU,r[rows]);
if vn~=v; break; endif;
v=vn;
end;
return vn;
endfunction
function xinv (A,n=20)
## Cumpute the invers of A using residual iteration.
## You can specify an epsilon eps with eps=... as last parameter.
if (isvar("eps")); return xlgs(A,id(cols(A)),eps=eps);
else (isvar("eps")); return xlgs(A,id(cols(A)));
endif
endfunction
function xlusolve (A,b,n=20)
## Compute the solution of Ax=b for LU-matrices A.
## E.g., A may be a lower triangle matrix.
## n is the number of residual iterations.
## The iteration stops, when there is no more change in x.
## You can specify an epsilon eps with eps=... as last parameter.
if (isvar("eps")); localepsilon(eps); endif;
v=lusolve(A,b);
loop 1 to n;
r=residuum(A,v,b);
vn=v-lusolve(A,r);
if vn~=v; break; endif;
v=vn;
end;
return vn;
endfunction
function xpolyval (p,t,n=20)
## Evaluate the polynomial at values t using exact residual
## operations. The iteration stops, when there is no more
## change in the value.
## You can specify an epsilon eps with eps=... as last parameter.
if (isvar("eps")); localepsilon(eps); endif;
s=t; si=size(t);
A=id(cols(p));
b=flipx(p)';
loop 1 to prod(si);
A=setdiag(A,-1,-t{#});
v=lusolve(A,b);
loop 1 to n;
r=residuum(A,v,b);
vn=v-lusolve(A,r);
if vn~=v; break; endif;
v=vn;
end;
s{#}=v[cols(p)];
end;
return s;
endfunction
function xeigenvalue (a,l,x=0)
## xeigenvalue(A,l) returns an improved the eigenvalue of A.
## l must be closed to a simple eigenvalue.
## Returns the eigenvalue and the eigenvector.
l1=l;
if x==0; x=kernel(a-l1*id(cols(a))); endif;
if x==0; error("Eigenvalue generation failed!"); endif;
x=x[:,1];
repeat;
x=(a-l1*id(cols(a)))\x;
x=x/sqrt(x'.x);
l2=residuum(x',residuum(a,x,0),0);
if l1~=l2; return {l2,x}; endif;
l1=l2;
end;
endfunction
|