Anonymous user
Gaussian elimination: Difference between revisions
m
→{{header|REXX}}: added/changed comments and whitespace, aligned statements better.
m (→{{header|REXX}}: added/changed comments and whitespace, aligned statements better.) |
|||
Line 4,034:
Only '''8''' (fractional) decimal digits were used for the output display.
<lang rexx>/*REXX program solves Ax=b with Gaussian elimination and backwards substitution. */
parse arg iFID . /*obtain optional argument from the CL.*/▼
numeric digits 1000 /*heavy─duty decimal digits precision. */
▲parse arg iFID . /*obtain optional argument from the CL.*/
if iFID=='' | iFID=="," then iFID= 'GAUSS_E.DAT' /*Not specified? Then use the default.*/
do rec=1 while lines(iFID) \== 0 /*read the equation sets.
#=
do $=1 while lines(iFID) \== 0 /*process the equation. */
z= linein(iFID);
if $==1 then do; say; say center(' equations ', 75, "▓"); say
end /* [↑] if 1st equation, then show hdr.*/
say z /*display an equation to the terminal. */
if left(space(z), 1)=='*' then iterate /*Is this a comment? Then ignore it.*/
#= # + 1;
do e=1 for n; a.#.e= word(z, e)
end /*e*/ /* [↑] process A numbers. */
b.#= word(z, n + 1)
end /*$*/
if #\==0 then call Gauss_elim /*Not zero? Then display the results. */
Line 4,054:
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
Gauss_elim:
do i=jp to n; _= a.j.j / a.i.j
do k=jp to n; a.i.k= a.j.k - _ * a.i.k
end /*k*/
b.i= b.j - _ * b.i
end /*i*/
end /*j*/
x.n= b.n / a.n.n
do j=n-1 to 1 by -1; _= 0
do i=j+1 to n; _= _ + a.j.i * x.i
end /*i*/
x.j= (b.j - _) / a.j.j
end /*j*/ /* [↑] uses backwards substitution. */
▲ numeric digits 8 /*for the display, only use 8 digits. */
say center('solution', 75, "═"); say /*a title line for articulated output. */
do o=1 for n; say right('x['o"] = ", 38) left('', x.o>=0) x.o/1
Line 4,137 ⟶ 4,136:
parse arg iFID . /*obtain optional argument from the CL.*/
if iFID=='' | iFID=="," then iFID= 'GAUSS_E.DAT' /*Not specified? Then use the default.*/
pad= left('', 23)
do rec=1 while lines(iFID) \== 0 /*read the equation sets. */
#=0 /*the number of equations (so far). */
do $=1 while lines(iFID) \== 0 /*process the equation. */
z= linein(iFID);
if $==1 then do; say; say center(' equations ', 75, "▓"); say
end /* [↑] if 1st equation, then show hdr.*/
say z /*display an equation to the terminal. */
if left(space(z), 1)=='*' then iterate /*Is this a comment? Then ignore it.*/
#= # + 1;
do e=1 for n; a.#.e= word(z, e); oa.#.e= a.#.e
end /*e*/ /* [↑] process A numbers; save orig.*/
b.#= word(z, n
end /*$*/
if #\==0 then call Gauss_elim /*Not zero? Then display the results. */
say
do i=1 for n; r=0 /*display the residuals to the terminal*/
do j=1 for n; r=r + oa.i.j * x.j /*
end /*j*/ /*
r= format(r
say pad 'residual['right(i, length(n) )"] = " left('', r>=0) r /*right justify.*/
end /*i*/
Line 4,162 ⟶ 4,161:
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
Gauss_elim:
do i=jp to n; _= a.j.j / a.i.j
do k=jp to n; a.i.k= a.j.k - _ * a.i.k
end /*k*/
b.i= b.j - _ * b.i
end /*i*/
end /*j*/
x.n= b.n / a.n.n
do j=n-1 to 1 by -1; _= 0
do i=j+1 to n; _= _ + a.j.i * x.i
end /*i*/
x.j= (b.j - _) / a.j.j
end /*j*/ /* [↑] uses backwards substitution. */
▲ numeric digits 8 /*for the display, only use 8 digits. */
say center('solution', 75, "═"); say /*a title line for articulated output. */
do o=1 for n; say right('x['o"] = ", 38) left('', x.o>=0) x.o/1
| |||