Help18:=proc(): print(`HWgE(K1,K2)`):end:
read `DMB.txt`:

#HWgE(K1,K2): tries out K2 different random parameters for HWg (the generalized Hardy-Weinberg)
#where thee entries of the matrix are from 1 to K1
#and outputs where the genotype gets extinct
#that give you a stable steady-state
HWgE:=proc(K1,K2)  local M,T,co,i,L:
co:=0:

for i from 1 to K2 do
M:=RandMat(3,K1):
T:=HWg(u,v,M);
L:=SSSgN(T,[u,v]);
if abs(L[1]-1)<0.001  or abs(L[2]-1)<0.001 or abs(L[1]-L[2])<0.001 then
co:= co+1:
fi:
od:
evalf(co/K2):
end: