Gale-Shapley Algorithm for R

Dear R-helpers,

I'm not a speciallist in writing complex functions, and the function still
very rusty (any kind of suggestions are very welcome). I want to implement
Gale-Shapley algorithm for R Language. It is based on 
http://www.jstor.org/stable/10.2307/2312726 Gale and Shapley (1962) , and
it has evolved to 
http://kuznets.fas.harvard.edu/~aroth/papers/Gale%20and%20Shapley.revised.IJGT.pdf
several applications  in many languages, C++, JAVA, 
http://stackoverflow.com/questions/2526042/how-can-i-implement-the-gale-shapley-stable-marriage-algorithm-in-perl
Perl , SQL, and so on. I manage to edit one version of it  to R.2.13.1.
So, I ask if it was allready implemented (I couldn't find any on R topic),
and if there is models and manners to make it more efficiently, add errors
check, options, etc.

At Berkeley's http://mathsite.math.berkeley.edu/smp/smp.html  MathSite 
there is a very straighfoward example of the algorithm and its steps.

My implementation follow the principle:

1. All men (or women) seeks for their best partner.

2. If there is no tie in a column (or row), stop.

3. If there is a tie, removes the worst-partners-tied and seek again the
second-best (till n-best) alternative.

The function is working right up to 6x6 matrices. But it needs a lot of
improvement.

Here it is the "gsa" function:

gsa(m, n, preference.row, preference.col, first)

###
Where:

m: for number of rows
n:  number of columns
preference.row: matrix with preference ordered in its positions by row
(see example).
preference.col: matrix with preference ordered in its positions by column
(see example).
first: Who is the first to propose (1 to men, 2 to women).
########

gsa <- function(m, n, preference.row, preference.col, first)
{
# m: number of rows (men)
# n: number of columns (women)
# first 1 for row (men); and 2 for column (women)
#
# Two Auxiliary functions:
# 1:
min.n <- function(x,n,value=TRUE){
s <- sort(x, index.return=TRUE)
if(value==TRUE){s$x[n]} else{s$ix[n]}}

# 2:

max.n <- function(x,n,value=TRUE){
s <- sort(x, decreasing=TRUE, index.return=TRUE)
if(value==TRUE){s$x[n]} else{s$ix[n]}}
#############################################################

s <- NULL
test_s <-NULL
loop <- 2 # O loop ? necess?rio a partir do 2.
step.1 <- matrix(0,ncol=n, nrow=m)
step.2 <- matrix(0,ncol=n, nrow=m)
store <- NULL
r <- NULL

# Men proposing first:

if (first==1)
        {
step.1 <- matrix(0,ncol=n, nrow=m)
for (i in 1:n)
                {
step.1[i,][preference.row[i,]==min.n(preference.row[i,],n=1)] <- 1
                }
for (i in 1:n){s[i] <- sum(step.1[,i])}
test_s <- s>1
while (any(test_s==TRUE)==TRUE)
                        {
if (any(test_s==TRUE)==TRUE) {
position1 <- which(s>1)
position2 <- vector('list')
position3 <- vector('list')
position4 <- NULL
position5 <- 1:m
for (k in 1:length(position1)){position2[[k]] <-
which(step.1[,position1[k]]==1)
position3[[k]] <-
which(preference.col[,position1[k]]>min(preference.col[position2[[k]],position1[k]]))
x <- which(position3[[k]]%in%position2[[k]])
position3[[k]] <- position3[[k]][x]
step.1[position3[[k]],position1[k]] <- 0}
for (t in 1:n){position4[t] <-
if(sum(step.1[,t])==0){0}else{which(step.1[,t]==1)}}
position4 <- position4[position4 >0]
position5 <- position5[-position4]
store <- append(position5, store)
r <- rle(sort(store))
for (j in
position5){step.1[j,][preference.row[j,]==r$lengths[r$values==j]+1] <- 1}
for (i in 1:n){s[i] <- sum(step.1[,i])}
test_s <- s>1
                                        }else{
step.1 <- matrix(0,ncol=m, nrow=n)
for (i in 1:m){step.1[i,][preference.row[i,]==min(preference.row[i,])] <-
1}
return(step.1)}
loop <- loop + 1
                        } #end of while
        }

# Women proposing first:

if (first==2)
        {
step.2 <- matrix(0,ncol=n, nrow=m)
for (i in 1:n)
                {
step.2[,i][preference.col[,i]==min.n(preference.col[,i],n=1)] <- 1
                }
for (i in 1:n){s[i] <- sum(step.2[i,])}
test_s <- s>1
while (any(test_s==TRUE)==TRUE)
                        {
if (any(test_s==TRUE)==TRUE) {
position1 <- which(s>1)
position2 <- vector('list')
position3 <- vector('list')
position4 <- NULL
position5 <- 1:m
for (k in 1:length(position1)){position2[[k]] <-
which(step.2[position1[k],]==1)
position3[[k]] <-
which(preference.row[position1[k],]>min(preference.row[position1[k],position2[[k]]]))
x <- which(position3[[k]]%in%position2[[k]])
position3[[k]] <- position3[[k]][x]
step.2[position1[k],position3[[k]]] <- 0}
for (t in 1:n){position4[t] <-
if(sum(step.2[t,])==0){0}else{which(step.2[t,]==1)}}
position4 <- position4[position4 >0]
position5 <- position5[-position4]
store <- append(position5, store)
r <- rle(store)
for (j in
position5){step.2[,j][preference.col[,j]==r$lengths[r$values==j]+1] <- 1}
for (i in 1:n){s[i] <- sum(step.2[i,])}
test_s <- s>1
                                        }else{
step.2 <- matrix(0,ncol=m, nrow=n)
for (i in 1:m){step.2[i,][preference.col[,i]==min(preference.col[,i])] <-
1}
step.2}
loop <- loop + 1
                        } # End of 2nd while
        }
if (first==1) {print(step.1)}
if (first==2) {print(step.2)}
}

#####################
# Here it goes one 4x4 example:

m <- seq(1:4)
n <- seq(1:4)
preference.row <- matrix(0,ncol=length(m), nrow=length(m))
preference.col <- matrix(0,ncol=length(n), nrow=length(n))

for (i in 1:length(m))
{
preference.row[i,] <- sample(m, size=length(m), rep=FALSE)
preference.col[,i] <- sample(n, size=length(n), rep=FALSE) # Note a
orienta??o por coluna!
}
gsa(m = 4, n = 4, preference.row = preference.row, preference.col =
preference.col, first=2)

# The result is a zero-one matrix which indicates blocking pairs.
############################################################

Thank you, and please let me know, any bugs and improvements.