Incorrect/inconsistent equality of custom units (value rounding error)

### Description
```
>>> import astropy.units as u
>>> Hz = u.def_unit('Hz', u.cycle / u.s)
>>> kHz = u.def_unit('kHz', 1000 * Hz)
>>> 1 * kHz == 1000 * Hz
True
>>> 1000 * Hz == 1 * kHz
False
```
Of course, both equalities should return True, and in general equalities should be commutative.

(As an aside, I understand why Astropy defines Hz as 1/s, but in my field of work, cycle/s is more convenient and less error-prone, so I want to be able to redefine it.)

### System Details
This was found on an Anaconda install on Windows:
- Windows-10-10.0.19042-SP0
- Numpy 1.20.3
- pyerfa 2.0.0
- astropy 5.0
- Scipy 1.7.1
- Matplotlib 3.5.0

Start by reproducing the custom-unit equality example through the astropy.units entry point, then trace how equality and value conversion handle the two operand orders. Done means both comparisons return True consistently, with regression coverage for the demonstrated Hz and kHz definitions.

Written by the indexing model from the issue text.