One nice way of doing this is using the new [LATERAL JOIN](https://www.postgresql.org/docs/current/static/queries-table-expressions.html) syntax in Postgres in conjunction with the <-> operator.
A lateral join essentially runs the right hand query for each row in the left hand one, and can avoid some quite nasty array logic and sub-queries, which you would otherwise need to use to find k nearest neighbours.

    CREATE TABLE sometable AS
    SELECT 
        poly.id, array_agg(pts.id) 
      FROM
         (SELECT id, geom FROM polygons) poly
      CROSS JOIN LATERAL
         (SELECT id, geom 
            FROM points pt
            ORDER BY pt.geom <-> poly.geom LIMIT 10
      ) pts;

There is a good [blog by Paul Ramsey](https://carto.com/blog/lateral-joins) explaining this technique.

There are a few gotchas though. The <-> operator operates on the spatial index, but generally requires that one of the two geometries is a constant. For small tables with LIMIT 1, this may not matter, but, as always EXPLAIN is your friend. I have found, in practice, that for very large tables (I have done this for tables of 50 million points, to find each points 6 nearest neighbours, which in the naive case gives you 50 x 50 million possible combinations) you might need to add  
    
     ST_DWithin(pt.geom, poly.geom, distance)

to the query, which will use the spatial index, but requires you to know beforehand a value for distance that is guaranteed to find at least 10 points for each polygon.

      SELECT 
        poly.id, array_agg(pts.id) 
      FROM
         (SELECT id, geom FROM polygons) poly
      CROSS JOIN LATERAL
         (SELECT id, geom 
            FROM points pt
            WHERE ST_DWithin(pt.geom, poly.geom, some_distance)
            ORDER BY pt.geom <-> poly.geom LIMIT 10
      ) pts;

I have found, in practice, that for fairly large tables where you are looking for several nearest neighbours, this can be done in a loop, where you start with a small search distaince in ST_DWithin, and gradually increase it for those polygons you haven't round the nearest 10 points for (which obviously requires a where clause in the initial select). Sadly, while it is syntactically very nice LATERAL JOINS and the <-> operator by themselves do not guarantee fast execution automatically, when you are looking for the k nearest neighbours for multiple geometries.

**Note**. It is possible to use **LEFT JOIN LATERAL** syntax too, but in this instance, you have to add **ON TRUE**, after the right-hand table, eg,

    SELECT 
      poly.id, array_agg(pts.id) 
    FROM
      (SELECT id, geom FROM polygons) poly
    LEFT JOIN LATERAL
      (SELECT id, geom 
        FROM points pt
        WHERE ST_DWithin(pt.geom, poly.geom, some_distance)
        ORDER BY pt.geom <-> poly.geom LIMIT 10
    ) pts ON TRUE;