Can someone share nicely written code for doing the same?

That is, given two coprime p,q find p.q^-1 mod 10e9+7.

From some time ago on Project Euler, I found this code to calculate x^{-1}

```
uint32_t invert( uint32_t x )
{
if (x==0) {
cerr << "Invert 0 not allowed";
exit(1);
}
if (x==1) return 1;
uint32_t a = invert(mod%x) * (uint64_t) (mod-mod/x) % mod;
return a;
}
```

The alternative is to rely on the fact that x^{p-1} \equiv 1, so that x^{-1} = x^{p-2} modulo p.

The alternative is true only when p is prime.