Demo Gauss-eliminatie voor Lineaire Algebra

Eerst moeten we een paar  hulpfuncties initialiseren (door Shift+Enter op deze groep)

generate[m_,n_,max_,pivots_]:=Module[{k,l,A,B,lp=Length[pivots]},
If[lp>m,Print["Number of pivots too big"];Return[Null]];
If[pivots[[lp]]>n-1,Print["pivotindex too big"];Return[Null]];
A=Table[RandomInteger[{-max,max}],{k,m},{l,n}];
For[k=1,k<=lp,k++,
For[l=1,l<pivots[[k]],l++,A[[k,l]]=0]
];
For[k=lp+1,k<=m,k++,
For[l=1,l<=n,l++,A[[k,l]]=0]
];
B=Table[0,{k,m},{l,m}];
While[Det[B]==0,B=Table[RandomInteger[{-4,4}],{k,m},{l,m}]];
Return[B.A];
];
gausstep[A_]:=Module[{m,n,k,l,B,test=0,vec},B=A;m=Length[B];n=Length[Transpose[B]];If[m==1,Return[B]];For[k=1,k≤n,k++,For[l=1,l≤m,l++,If[Abs[B[[l,k]]]≠0,Goto[gotyou]]]];Return[B];Label[gotyou];If[l>1,vec=B[[l]];B[[l]]=B[[1]];B[[1]]=vec;Return[B]];test=(∑_(l=2)^m)Abs[B[[l,k]]];If[test>0,For[l=2,l≤m,l++,B[[l]]=B[[l]]-((B[[l,k]] B[[1]])/(B[[1,k]]))];Return[B]];B=Join[{B[[1]]},gausstep[Drop[B,1]]];Return[B]];
jordanstep[A_]:=Module[{B,m,n,k,l,r},
B=A;
m=Length[B];
n=Length[Transpose[B]];
For[k=m,k≥1,k--,
For[l=1,l≤n,l++,
If[B[[k,l]]≠0,If[B[[k,l]]≠1,B[[k]]=((B[[k]])/(B[[k,l]]));Return[B]];
If[B[[k,l]]==1,If[(∑_r^(k-1))Abs[B[[r,l]]]≠0,For[r=1,r<k,r++,B[[r]]=B[[r]]-B[[r,l]] B[[k]]];Return[B],l=n+1];
];
];
];
];
Return[B];
]
MyMatrixForm[A_]:=Module[{B=Transpose[A],K,X,len,k},len=Length[B]-1;X=Transpose[{Table[x[k],{k,len}]}];K=Transpose[{Last[B]}];B=Transpose[Drop[B,-1]];B=Row[{MatrixForm[B],MatrixForm[X],"  =  ",MatrixForm[K]}];Return[B]]

Hier is een aardige manier om interessante stelsels te produceren

A=generate[5,7,10,{1,2,4}];

A[[1,1]]=0;

Clear[b];A[[5,7]]=b;

A//MyMatrixForm

Vervolgens kunnen we de Gauss-eliminatie doen door herhaald het volgende commando te geven

A=gausstep[A];A//MyMatrixForm

Tenslotte een Jordan reductie

A=jordanstep[A];A//MyMatrixForm