Calculating p-value using R i want to calculate the p-value of a statistics test using R. I'm aware of the existing function t.test(x) but what i want to do is to determine the p value through monte carlo simulations. Therefore i defined this code:
mean.test <- function(x, y, B=10000,
alternative=c("two.sided","less","greater"))
{

p.value <- 0
alternative <- match.arg(alternative)

s<-replicate(B, (mean(sample(c(x,y), B, replace=TRUE))-mean(sample(c(x,y), 
B, replace=TRUE)))) # random samples of test statistics
t <- mean(x) - mean(y) #teststatistics t    
p.value <- 2 * (1- pnorm(mean(s)))   #try to calculate p value 

data.name <- deparse(substitute(c(x,y)))
names(t) <- "difference in means"
zero <- 0
names(zero) <- "difference in means"
return(structure(list(statistic = t, p.value = p.value,
method = "mean test", data.name = data.name,
observed = c(x,y), alternative = alternative,
null.value = zero),
class = "htest"))
}

When running 
> set.seed(0)
> mean.test(rnorm(1000,3,2),rnorm(2000,4,3))

this is supposed to return
     mean test
data: c(rnorm(1000, 3, 2), rnorm(2000, 4, 3))
difference in means = -1.0967, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0

this return this though:
        mean test
 data:  c(rnorm(1000, 3, 2), rnorm(2000, 4, 3))
 difference in means = -1.0967, p-value = 0.9999
 alternative hypothesis: true difference in means is not equal to 0

What's the error ?
 A: Please see the caveats in my comments to the question.  This is just a conceptual example. It is not meant to be presented as a statistically sound test.
The following code will compare the means of two samples, and produce a p-value by resampling.  Note that if the samples are normally distributed and there are sufficient replications, that the p-value is very close to that of the t.test function.
The heart of the function is very simple.  I added some code to produce a nicer output.
Mean.test = function(X, Y, r=1000, digits=4){
  XY = c(X, Y)
  Diff = abs(mean(X) - mean(Y))
  Count = 0
  for(i in 1:r){
  S1 = sample(XY, size=length(X), replace=TRUE)
  S2 = sample(XY, size=length(Y), replace=TRUE)
  Diff.s = abs(mean(S1) - mean(S2))
  if(Diff.s >= Diff){Count=Count+1}
  }
  P.value = Count/r
  Z=data.frame(
    Statistic= rep("Boogida", 3),
    Mean    = rep(NA, 3),
    SD      = rep(NA, 3),
    n       = rep(NA, 3),
    stringsAsFactors=FALSE)

  Z[1,] = c("Group 1", signif(mean(X),digits),signif(sd(X),digits),
            signif(length(X),digits))
  Z[2,] = c("Group 2", signif(mean(Y), digits), signif(sd(Y), digits),
            signif(length(Y),digits))
  Z[3,] = c("Difference", signif(abs(mean(X) - mean(Y)), digits), NA, NA)

  colnames(Z)[1]=""

  U = data.frame(
      Statistic = rep("Boogida", 1),
      p.value    = rep(NA, 1),
      stringsAsFactors=FALSE)

  colnames(U)[1]=""

  U[1,] = c("p-value", signif(P.value, digits))

  W = list(Summary=Z, 
           Result=U)

  return(W)
}

A = rnorm(n=19, mean=6, sd=3)
B = rnorm(n=21, mean=8, sd=3)

Mean.test(A, B, r= 10000)

t.test(A, B)

