Ecco mostrato il testo dell'algoritmo probabilistico. I simboli sottolineati sono dei link.

n = value read "n? "

k = value read "k? "

N = binomial (n+1,k+1)

S = S = QQ [apply (toList((0,0)..(k,n-k-1)), i -> x_(toSequence i))]

X = vars S

M = X -> id_(S^(k+1))|transpose genericMatrix(S,X_0_0,n-k,k+1)

par = X -> exteriorPower(k+1,M(X))

G = matrix({apply(N,i->(flatten entries par(X))_i)})

J = jacobian G

f = map(S,QQ)

U = l -> ( a=l_0; for i from 1 to #l-1 do ( a=a||(l_i); ); a)

m = p ->substitute(M(X),matrix(S,{flatten entries p}))

par' = p -> exteriorPower(k+1,p)

s = value read "s? "

punti = apply(s,i->f(random(QQ^(k+1),QQ^(n-k))))

MP = apply(s,i->m(punti_i))

PAR = apply(s,i->par’(MP_i))

time P = apply(s,j->matrix(S,{apply(N,i->(flatten entries PAR_j)_i)}))

JPpar = apply(s,i->substitute(J,matrix{flatten entries punti_i}))

Ja = apply(s,i->(JPpar_i)||(P_i))

rank U(Ja).

L'output rank(U(Ja)) è la dimensione cercata della varietà s- secante alla Grassmanniana G(k,n).