Multinomial regression with restriction

Modeling
1

Dear Stan Users,

I am new in stan and I would like to fit a multinomial model with a restriction on the probabilities. I believe it is different from a softmax model, but I am not really sure about it.
Below, you can see my stan model and R code. However, when compiling in R a got some error messages. I looked for help on internet, but no success so far. Please, let me know if some of you already work with a similar model and how to code it in stan.

Thanks.

data {
  int K;
  int N;
  int D;
  int y[N, K];
  matrix[N, D] X;
}
parameters {
  real beta1[3];
  real beta2[3];
}
transformed parameters {	
  real linpred1[N];
  real linpred2[N];
  matrix[N, D] prob;
  for(i in 1:N) {
  linpred1[i] = exp(X[i,1]*beta1[1] + X[i,2]*beta1[2] + X[i,3]*beta1[3]);
  linpred2[i] = exp(X[i,1]*beta2[1] + X[i,2]*beta2[2] + X[i,3]*beta2[3]);
  prob[i,2] = linpred1[i]/(1+linpred1[i] + linpred2[i]);
  prob[i,3] = linpred2[i]/(1+linpred1[i] + linpred2[i]);
  prob[i,1] = 1 - prob[i,2] - prob[i,3];
  }
}
model {
    for (i in 1:N) {
      y[i,] ~ multinomial(to_row_vector(prob[i,]) );
    }
} 

set.seed(123)

# Generating the probabilities ---------------------------------------
x1 <- seq(-1, 1, l = 100)
x2 <- rbinom(100, prob = 0.5, size = 1)
X <- model.matrix(~ x1 + x2)
beta1 <- c(0, 0.8, -1)
beta2 <- c(0.5, -0.8, -1)

eta1 <- X%*%beta1
eta2 <- X%*%beta2

pi1 <- exp(eta1)/(1+ exp(eta1) + exp(eta2))
pi2 <- exp(eta2)/(1+ exp(eta1) + exp(eta2))
pi0 <- 1 - pi1 - pi2

# Simulating data set ------------------------------------------------
Y <- matrix(NA, ncol = 3, nrow = 100)
for(i in 1:100) { 
  Y[i,] <- rmultinom(1, size = 1, prob = cbind(pi0, pi1, pi2)[i,])
}

library(rstan)
data <- list("y" = Y, "K" = 3, "D" = 3, "N" = 100, "X" = X, 
             beta1 = c(0.1, 0.6, 0), beta2 = c(0.1, 0.8, 0))

fit1 <- stan(
  file = "mult1.stan",  # Stan program
  data = data,    # named list of data
  chains = 4,             # number of Markov chains
  warmup = 1000,          # number of warmup iterations per chain
  iter = 2000,            # total number of iterations per chain
  cores = 2,              # number of cores (could use one per chain)
  refresh = 0             # no progress shown
)