You can get a lot of the way (just not to `exp((E_old - E_new) / Temp)`) by explaining the origin of the algorithm:
People who work metals have a process for making metals more supple, called "annealing". In annealing, you heat the metal up almost to the point where it loses its shape, then let it cool slowly. In the atomic picture, this chunk of metal has knots of stress and tiny mis-aligned crystals, but if we let the atoms jiggle, the knots can just work themselves out and the crystals can grow, because atoms want to align with their neighbors. It works the other way too though: If atoms jitter randomly then at any time they can exhibit a chaotic state, so if you suddenly stop that jitter ("lowering the temperature too fast"), then you just get a random state. So there is a sweet spot. You want just enough jitter to undo knots, but not so much jitter that you create new ones. By slowly reducing the temperature, you get to pass some special points where the neighbors have about the same impact on the crystal as the jitter, and as you slowly pass that, the jitter starts undoing the knots for you and then growing the crystals for you. You still often don't get perfection, but you can get closer.
What does this have to do with optimization? Physicists like to think of this "aligning with neighbors" process as minimizing some "interaction energy". But you can minimize whatever function you want with the same approach: (1) take any state you want, (2) let it randomly jitter to neighboring states, (3) accept the new states when it makes your function smaller, but also (with some probability p) when it makes your function larger, and (4) slowly lower p as you allow the system to randomly jitter around.
You can also get a good idea of what's going on with a graph:
\ _b / _
\_/ \ / _|-- B
\_/ _|-- A
Let's explain simulated annealing by looking that the barrier above (b), which has two different energy costs to jumping over it: B from the left and A + B from the right.
Here we understand the jitter as having a typical energy E which lets it jump up the graph, but it still tends to fall in the holes. When E > A + B then it will jump from the deeper valley on the right into the left and back. But as we get into the range B < E < A + B, then it's possible to jump from left to right, but not easily from right to left. The more we allow the system to jitter in this region as we slowly reduce E, the more likely the state ends up on the right somewhere, and then from that point on the algorithm mirrors your normal gradient descent algorithm.
> You can get a lot of the way (just not to `exp((E_old - E_new) / Temp)`) by explaining the origin of the algorithm
The exp is from the so-called Boltzmann factor [1] which (in somewhat simplified terms) describes the likelihood that a physical system in equilibrium at temperature T will transition from a state with energy E_old to one with energy E_new. So sure, you can derive the whole algorithm from classic "annealing" and physics alone.
I mean, I know where it comes from; I just don't think that there's a nice hand-wavy way to get people to understand what's going on there: you really need to start with how uncertainties multiply across a system, define entropy as log(W) in the NVE ensemble, discover temperature hiding there too, then switch to the NVT ensemble (i.e. connect to an infinite reservoir of energy at a certain temperature) -- then you can get at Boltzmann factors.
Even then, as was pointed out to me recently, it's not clear that the "always accept a lower energy, accept a higher energy proportional to the Boltzmann factor" yields a steady-state Boltzmann distribution -- especially because that would imply that the result is independent of the number of configurations with a given energy (the "density of states"), which seems surprising. If that doesn't hold then it's either a convenient approximation or a bit of magical thinking...
People who work metals have a process for making metals more supple, called "annealing". In annealing, you heat the metal up almost to the point where it loses its shape, then let it cool slowly. In the atomic picture, this chunk of metal has knots of stress and tiny mis-aligned crystals, but if we let the atoms jiggle, the knots can just work themselves out and the crystals can grow, because atoms want to align with their neighbors. It works the other way too though: If atoms jitter randomly then at any time they can exhibit a chaotic state, so if you suddenly stop that jitter ("lowering the temperature too fast"), then you just get a random state. So there is a sweet spot. You want just enough jitter to undo knots, but not so much jitter that you create new ones. By slowly reducing the temperature, you get to pass some special points where the neighbors have about the same impact on the crystal as the jitter, and as you slowly pass that, the jitter starts undoing the knots for you and then growing the crystals for you. You still often don't get perfection, but you can get closer.
What does this have to do with optimization? Physicists like to think of this "aligning with neighbors" process as minimizing some "interaction energy". But you can minimize whatever function you want with the same approach: (1) take any state you want, (2) let it randomly jitter to neighboring states, (3) accept the new states when it makes your function smaller, but also (with some probability p) when it makes your function larger, and (4) slowly lower p as you allow the system to randomly jitter around.
You can also get a good idea of what's going on with a graph:
Let's explain simulated annealing by looking that the barrier above (b), which has two different energy costs to jumping over it: B from the left and A + B from the right.Here we understand the jitter as having a typical energy E which lets it jump up the graph, but it still tends to fall in the holes. When E > A + B then it will jump from the deeper valley on the right into the left and back. But as we get into the range B < E < A + B, then it's possible to jump from left to right, but not easily from right to left. The more we allow the system to jitter in this region as we slowly reduce E, the more likely the state ends up on the right somewhere, and then from that point on the algorithm mirrors your normal gradient descent algorithm.