To expand on Kevin's comment (and using an answer since a comment doesn't have enough characters!) : one other obvious-but-relevant constraint that's going to be an issue for large values of $N$ is that the number of possible answers that can be acheived is a function strictly of $n$ and not of the $a_i$ (modulo a few minor issues like repetition of values etc): there are only $C_{n-1}$ 'operation' trees with $n$ leaves, where $C_n$ are the Catalan numbers, approximately $4^n$ plus some polynomial factors; since each of the $(n-1)$ internal nodes can be filled with one of 4 binary operators that adds another factor of $4^{n-1}$ to the total; and of course the $a_n$ can be permuted in $n!$ ways, so the overall bound is something like $n^n\alpha^n$ for $\alpha \approx 16/e \approx 5.9$; concretely, there are a maximum of $C_3 * 4^3 * 4!$ = 5 * 64 * 24 = 7680 possibilities for the $n=4$ case, so if $N$ is larger than this, you're guaranteed to have gaps regardless of the values for $a_i$ (and as $N$ goes to infinity, the number of possible solutions won't change, so the probability will be approximately $c/N$ for some constant $c$ - one interesting question would be what $c$ is likely to be, or more specifically, how many 'collisions' will reduce the number below the hard bound? There are a lot of obvious cases (e.g., swapping the terms of a summand won't make any difference) so the number of possible values is well below the explicit upper bound, but other than wildly speculating that it's still superexponential in $n$ I'm not even sure how to begin calculating it.

One conjecture that seems both plausible and accessible to me is that any 'numeric' bound is likely to hew closely to the 'algebraic' upper bound (where you treat the $a_i$ as independent variables rather than numbers and consider only the equivalence of the resulting functions and not the values they evaluate to); more specifically, if $F_n$ is the aforementioned number of algebraically distinct functions for a given $n$, then $\Sigma_{a_0=0}^{N}\Sigma_{a_1=0}^{N}\dots\Sigma_{a_n=0}^{N}(\sharp values)= F_nN^n(1-o(1))$, where $(\sharp values)$ is the number of distinct values acheivable for that set of $a_i$.