#c10.txt
Help10:=proc():print(` RT(n), IsStable1(M,W,pi,i1,i2), IsStable(M,W,pi), FindStable(M,W,K) `):end:
with(combinat):

#RT(n): Generates a pair of random preferance M,W, two lists-of-lists
#M[i][j]:= The ranking of Ms. j in the eyes of Mr. i
#W[j][i]:= The ranking of Mr. i  in the eyes of Ms. j

RT:=proc(n) local i,j:
[seq(randperm(n),i=1..n)],[seq(randperm(n),j=1..n)]:
end:

# i1 is married to j1
#i2 is married to j2
#i1 would rather be married to j2 and j2 would rather be married to i1
IsStable1:=proc(M,W,pi,i1,i2) local j1,j2:
j1:=pi[i1]:
j2:=pi[i2]:

  not(M[i1][j2]<M[i1][j1] and W[j2][i1]<W[j2][i2]) :
  
end:


#IsStable(M,W,pi): inputs preferenace lists-of-lists M,W, and a proposed matching is it stable?
#If it is not it outputs false, and an unstable pair i1,i2:
#
IsStable:=proc(M,W,pi) local n,i1,i2:

n:=nops(pi):

for i1 from 1 to n do
 for i2 from 1 to n do
  if i1<>i2 then
   if not IsStable1(M,W,pi,i1,i2) then
      RETURN(false, [i1,i2]):
   fi:
 fi:
od:
od:
true,true:
end:

#FindStable(M,W,K): Given the Men's and Women's preferance tables, starts with a random matching
#and keeps "fixing it" until hopefully you get a stable matching, we limit the number of
#iterations to K, if none found it will return FAIL. Otherwise it will return a stable
#mathing, followed by the number if iterations it took
FindStable:=proc(M,W,K) local n,pi,co,a,i1,i2,j1,j2:
print(`Warning: This is a stupid program that usually FAILS, for larger n`):
n:=nops(M):
pi:=randperm(n):

for co from 1 to K do
 a:=IsStable(M,W,pi):
 if a[1] then
   RETURN(pi,co):
 else
 i1:=a[2][1]:
 i2:=a[2][2]:
  j1:=pi[i1]:
  j2:=pi[i2]:
 if i1<i2 then
  pi:=[op(1..i1-1,pi),j2,op(i1+1..i2-1,pi),j1,op(i2+1..n,pi)]:
else
pi:=[op(1..i2-1,pi),j1,op(i2+1..i1-1,pi),j2,op(i1+1..n,pi)]:
fi:
fi:
od:
FAIL:
end: