"sterile.male"<-
function(a = 0.00, lambda = 4, mu = 1, N = 100, initial.size = 20,sterile=100,K=100)
{
	#
	# To set the random number seed, un-comment the following line 
	# and plug in any non-negative integer in [0,1023].
	# set.seed(123)
	#
	# 
	S <- rep(0, N)
	W <- S
	X <- S
	X[1] <- initial.size
	W[1] <- S[1]
	#
	# Note that this is time = 0.
	#
	for(i in 2:N) {
		if(X[i - 1] == 0) {
			if(a > 0) {
				# population is saved by immigration
				S[i] <- rexp(1, rate = a)
				W[i] <- W[i - 1] + 1
				X[i] <- X[i - 1] + 1
			}
			else {
				S[i] <- 0
				W[i] <- W[i - 1]
				X[i] <- X[i - 1]
			}
		}
		else {
			if(X[i-1]<K){
				#
				# Simulate a sojourn (inter-event) time, and compute the waiting time.
				#
				overall <- a + lambda * X[i - 1]*X[i-1]/(X[i-1]+sterile) + mu * X[i - 1]
				S[i] <- rexp(1, rate = overall)
				W[i] <- W[i - 1] + S[i]
				u <- runif(1)
				if(u < (mu * X[i - 1])/overall) {
					#event is a death
					X[i] <- X[i - 1] - 1
					}
				else {
					#event is a birth
					X[i] <- X[i - 1] + 1
				}
			}#end if X[i-1]<K
			else{#population at carrying capacity: event is a death
				S[i] <- rexp(1,rate=mu*X[i-1])
				W[i] <- W[i-1]+S[i]
				X[i] <- X[i-1]-1
			}
			}# end else
	}
	# end for loop on i
	#
	#
	plot(W,  X, type = "n", xlab = "Time", ylab = 
		"Population Size")
	lines(W, X, type = "s")
	return()
}