Consider the code for generating gamma distribution using MCMC in R:

N <- 2999
library(MCMCpack)

my_gamma <- function(x){
    return(x**N*exp(-2*x))
}

samples = MCMCmetrop1R(fun=my_gamma, theta.init=1,V=as.matrix(1))
print(samples)

The code says that the function returns Inf, although the mean is 1500, clearly in feasible range. I suspect that this might be due to large number computation in R.
For example, print(1500**3000) yields Inf. Is there a way to modify such computation limit to make the above code work?