Allow me to clarify a few points that I think have been confused in this thread.

In Stan’s implementation of dynamics Hamiltonian Monte Carlo a numerical trajectory is expanded in stages, randomly choosing a time direction and expanding the previous trajectory in that direction to double the number of states until a termination criterion is triggered. How that trajectory is used depends on what termination criterion is triggered.

(1) If a divergence (rapid change in density) is encountered in any of the new states then the states in the trajectory expansion are ignored and a sample is drawn from the previous trajectory.

(2) If the generalized No-U-Turn criterion is violated anywhere amongst the new states then the states in the trajectory expansion are ignored and a sample is drawn from the previous trajectory.

(3) If the generalized No-U-Turn criterion is violated for the entire expanded trajectory then a sample is drawn from the new expanded trajectory.

This procedure ensures the correct invariant distribution so long as the target density function is consistent, which is guaranteed by the confines of the Stan Modeling Language (up to errors in floating point arithmetic and the implementation of transcendental functions, etc). This then always ensures at least asymptotic consistency.

In the Stan Modeling Language `reject`

is implemented by throwing a domain error in the density function/gradient calculation which Stan’s Hamiltonian Monte Carlo implementation reinterprets as an infinite density. In other words if `reject`

is called then the Stan program will return positive infinity.

The jump from any finite target density to an infinite target density then triggers the divergence condition, which then triggers condition (1) above. All of the states in the current expansion, including the one that triggered the rejection, are thrown away and the sampler chooses a new state from previous trajectory.

If the `reject`

statement is within a conditional statement then it may not be called until the trajectory has already grown a bit, in which case the sampler will not necessarily return the initial point over and over again. If the `reject`

statement is always called, however, then the first expansion away from the initial point should be immediately rejected and sampler will return the same initial point over and over again.

Importantly Stan’s Hamiltonian Monte Carlo algorithm is *not* a Metropolis-Hastings algorithm and so any intuition and expectations based on accepting/rejecting a single proposed state won’t make much sense.