Have you measured the running time of this method? One concern is that the for loop would make the backprop time grow linearly with the dimension of the output. I’m wondering whether there is way to directly get the Jacobian without the for loop.