Quantum trajectories

Here’s your problem: You’re referencing integrator.u, not the u that’s provided separately as an argument. This way you fail to ever trigger the callback, so your trajectory doesn’t have any jumps. You should define the condition as follows:

function jump_condition(u, t, integrator)
    return sum(abs2, u) - integrator.p.threshold
end

I highly recommend plotting as a debugging strategy here. It would make the lack of jumps immediately obvious and point you toward the location of the problem. In general, it’s good practice to do some sanity checks and plots of a single trajectory before moving on to ensembles.

Btw., the line prob.p.threshold = rand() in prob_func is redundant. The remake call is sufficient.

(Further explanation of the issue: when the callback mechanism sees a sign flip in the condition, it uses interpolation to locate the exact crossing. The interpolated state values are not saved to integrator.u, but only provided through the argument u, so a state-dependent callback must use the latter.)