In the world of chess engines, what you're looking for is called multi-PV (multiple principal variations): telling the engine to look for and return up to n best moves and the associated principal variations for each. ("Up to" because sometimes there's only 1 legal move in a position.)
In python-chess, this is part of the engine protocols, as one of the optional parameters when calling
engine.analyse(). From the documentation:
coroutine analyse( [...] )
- multipv – Optional. Analyse multiple root moves. Will return a list of at most multipv dictionaries rather than just a single info dictionary.
.analyse() method tells the engine to analyse the position, returning a dictionary of information, which includes among other things the moves of the principal variation. Specifying for instance
.analyse(multipv = 3) will tell the engine to return a list of (up to) 3 such dictionaries -- the first moves in the
pv value of each dictionary should be what you are looking for.