Si vous pouvez limiter la distance maximale entre vos villes et votre position locale, profitez du fait qu'une minute de latitude (nord - sud) équivaut à un mile nautique.
Mettez un index sur votre table de latitude.
Créez-vous une fonction stockée haversine(lat1, lat2, long1, long2, unit) à partir de la formule haversine indiquée dans votre question. Voir ci-dessous
Ensuite, faites ceci, étant donné mylatitude, mylongitude et mykm.
SELECT *
from cities a
where :mylatitude >= a.latitude - :mykm/111.12
and :mylatitude <= a.latitude + :mykm/111.12
and haversine(:mylatitude,a.latitude,:mylongitude,a.longitude, 'KM') <= :mykm
order by haversine(:mylatitude,a.latitude,:mylongitude,a.longitude, 'KM')
Cela utilisera une boîte englobante de latitude pour exclure grossièrement les villes trop éloignées de votre point. Votre SGBD utilisera une analyse de plage d'index sur votre index de latitude pour sélectionner rapidement les lignes de votre table de villes qui valent la peine d'être prises en compte. Ensuite, il exécutera votre fonction haversine, celle avec tous les calculs sinus et cosinus, uniquement sur ces lignes.
Je suggère la latitude car la distance au sol de la longitude varie avec la latitude.
Notez que c'est grossier. C'est bien pour un chercheur de magasin, mais ne l'utilisez pas si vous êtes un ingénieur civil - la terre a une forme elliptique et cela suppose qu'elle est circulaire.
(Désolé pour le nombre magique 111.12. C'est le nombre de km dans un degré de latitude, c'est-à-dire en soixante milles nautiques.)
Voir ici pour une fonction de distance utilisable.