model <- function(t, state, parms) {
  with(as.list(c(state,parms)), {
    R <- s/(d + c1*N1 + c2*N2)
    dN1 <- N1*(b1*c1*R/(h1+c1*R) - d1)
    dN2 <- N2*(b2*c2*R/(h2+c2*R) - d2)
    return(list(c(dN1, dN2)))  
  }) 
}

# N1 is r-selected and N2 is K-selected:
s <- c(N1=0.1,N2=0.1)
p <- c(b1=1,b2=0.2,d1=0.5,d2=0.02,c1=2,c2=0.5,h1=0.1,h2=0.3,s=1,d=1)
plane(xmax=15,ymax=40,portrait=TRUE)

Rstar1 <- with(as.list(p),h1/(c1*(b1/d1-1)))
Rstar2 <- with(as.list(p),h2/(c2*(b2/d2-1)))
R01 <- with(as.list(p),b1/d1)
R02 <- with(as.list(p),b2/d2)
hR01 <- with(as.list(p),(b1/d1)*(c1*s/d)/(h1+c1*s/d))
hR02 <- with(as.list(p),(b2/d2)*(c2*s/d)/(h2+c2*s/d))
print(c(Rstar1=Rstar1,Rstar2=Rstar2,R01=R01,R02=R02,hatR01=hR01,hatR02=hR02))
# r-selected wins despite lower R0 and K