pytorch absolutely does not use the dual number formulation - there are absolutely no magic epsilons anywhere in pytorch's (or tensorflow's) code base. what you're calling duals are the adjoints where are indeed stored/cached on every node in pytorch graphs.
there's a reason no one uses dual numbers (non-standard analysis) for anything (neither autodiff nor calculus itself): because manipulating infinitesmals like this is fraught formal manipulation (it's algebra...) where as limits are much more rigorous (bounds, inequalities, convergence, etc.). my favorite question to ask the non-standard analysis n00bs is: please tell me under what conditions this is true
(dx/dy)(dy/dz)(dz/dx) = 1
edit:
anyone that thinks i'm wrong and this other guy is right should go and do some reading, eg where this guy tried to make this same point and got shot down:
https://math.stackexchange.com/a/341550
spoiler alert: there's a reason you had to learn epsilon-delta proofs and limits and it's not because your math professors are mean.
this is why i hate this kind of "TIL, gee whiz" math tidbits - they're full of exclamation marks and fancy sounding words ("non-archimedean rings" oooo fancy) but almost always come from a wikipedia level understanding, not actual research.