If your constraint is:
a \le D\beta \le w
Then if you multiply all by the inverse of D:
D^{-1}a \le D^{-1}D\beta \le D^{-1}w
This simplifies to:
D^{-1}a \le \beta \le D^{-1}w
Stan has overloaded the division operator /
for matrices/vectors to perform this multiplication by an inverse in more numerically stable and efficient way. So you would express your model in Stan as:
data {
int<lower=0> N; // Number of observations
int<lower=0> K; // number of explanatory variables.
matrix[N,K] X; // Matrix of explanatory variables
matrix[K,K] D; // Matrix that represents inequality constraints
vector[N] y; // Vector of dependent variables
vector[K] a; // Vector of lower inequality constraints
vector[K] w; // Vector of upper inequality constraints
}
transformed data {
vector[K] D_div_a = D / a;
vector[K] D_div_w = D / w;
}
parameters {
vector<lower=D_div_a, upper=D_div_w>[K] b;
real<lower=0> sigma;
}
model {
// More efficient likelihood for linear regression
y ~ normal_id_glm(X, 0, b, sigma);
}