May and Coron's result is a deterministic reduction from knowing a multiple of $\varphi(n)$, in this case $de-1$, to the factorization of $n$. This method involves lattice reduction and Coppersmith's algorithm.
But there is a much simpler probabilistic method to factor $n$ discovered by Miller, of the famed Miller-Rabin primality test. Suppose you are able to find a nontrivial square root of $1$, i.e.,
$$
x^2 \equiv 1 \pmod{n}, \quad x \neq \pm 1\,.
$$
Then $x^2 - 1 \equiv (x+1)(x-1) \pmod{n}$, and we find a factor by computing $\gcd(x \pm 1, n).$ But how do we find such a root? By exploiting the relation
$$
a^{\varphi(n)} \equiv 1 \pmod{n}\,.
$$
Let $t = \varphi(n) / \gcd\left(\varphi(n), 2^{\log_2(n)}\right)$ be $\varphi(n)$ with every power of 2 removed from it.
If we knew $\varphi(n)$ exactly, all we had to do was compute $a^{t} \bmod n$ for a random base $a$, and repeatedly square until we obtained a $1$—that value would be our $x$ as above. This will yield a factor depending on the order of $a$ modulo $p-1$ and $q-1$. This process naturally extends to any multiple of $\varphi(n)$, and in practice only a couple of bases are necessary to factor $n$.
Here's some Sage code that implements this:
def factor_from_d(n, d, e):
kphi = d*e - 1
# remove powers of 2 from phi multiple
kphi = kphi // gcd(kphi, 2^int(log(n, 2)))
while True:
# random base
b = randint(0, n)
x = power_mod(b, kphi, n)
# try to find a nontrivial square root of 1
while x != 1 and x != n - 1:
# found one!
# x^2 = 1 (mod n)
# x^2 - 1 = (x + 1)*(x - 1) (mod n)
if x^2 % n == 1:
return gcd(n, x + 1), gcd(n, x - 1)
x = x^2 % n
p = random_prime(2^512)
q = random_prime(2^512)
n = p*q
e = 2^16 + 1
d = inverse_mod(e, (p-1)*(q-1))
factor_from_d(n, d, e)
Coming back to May and Coron, their algorithm is conceptually simple, but more complex to implement from scratch. Here we want to find a small root of
$$
f(x) = n - x \pmod{de - 1}\,,
$$
which is guaranteed to exist because $n$ and $\varphi(n)$ are relatively close together. This root is precisely $p+q-1 = pq - (p-1)(q-1)$. Having $p+q$, we can solve the polynomial $x^2 - (p+q)x + n$, which has $p$ and $q$ as its roots.
Using Sage's small_roots function, we can do this easily:
p = random_prime(2^512)
q = random_prime(2^512)
n = p*q
e = 2^16 + 1
d = inverse_mod(e, (p-1)*(q-1))
P.<X> = Zmod(d*e-1)[]
f = X - n
pq = ZZ(f.small_roots(X=2^513, beta=0.9, epsilon=0.1)[0] + 1)
disc = pq^2 - 4*n
(pq + isqrt(disc))/2, (pq - isqrt(disc))/2
