Postgis Functional Index on st_pointonsurface() Ignored

Question

If I define the following table...

CREATE TABLE tbl_vast_importance(
  id bigserial PRIMARY KEY,
  geom geometry(MultiPolygon,4326)
);
CREATE INDEX ON tbl_vast_importance USING GIST (geom);
CREATE INDEX ON tbl_vast_importance USING GIST (st_pointonsurface(geom));

...and run the following query

SELECT

  id,
  geom

FROM

  tbl_vast_importance a
  JOIN tbl_vast_importance b ON (st_intersects(st_pointonsurface(a.geom),b.geom))

WHERE

  a.id != b.id

...the spatial index is not used; yet using st_intersects() on the geometry directly results in the geometry index being used properly. How do I force Postgres to use the proper index for st_pointonsurface()? It would be wasteful to create an explicit centroid column for this purpose.

Explain Results:

Nested Loop  (cost=0.00..1034058179737.31 rows=609220605748 width=784)
  Join Filter: ((a.id <> b.id) AND st_intersects(st_pointonsurface(a.geom), b.geom))
  ->  Seq Scan on tbl_vast_importance a  (cost=0.00..78073.11 rows=1351911 width=392)
  ->  Materialize  (cost=0.00..153484.67 rows=1351911 width=392)
        ->  Seq Scan on tbl_vast_importance b  (cost=0.00..78073.11 rows=1351911 width=392)
Planning time: 0.188 ms

postgis_full_version = "POSTGIS="2.2.0dev r12887" GEOS="3.5.0dev-CAPI-1.9.0 r3995" PROJ="Rel. 4.8.0, 6 March 2012" GDAL="GDAL 1.11.0, released 2014/04/16" LIBXML="2.9.1" LIBJSON="0.11.99" TOPOLOGY RASTER"

Answer 1

The reason that ST_Intersects(a.geom, b.geom) uses the index while ST_Intersects(ST_PointOnsurface(a.geom),b.geom)) does not is because while ST_Intersects will always use an index, if available, performing a bounding box test first, ie, &&, there is no way to know what the result of ST_PointOnSurface will be, so the optimizer no longer regards this as something to automatically apply an index to, based on pre-existing geometry bounding boxes -- as your EXPLAIN shows. The fact that your index uses the same function name is irrelevant -- one is a functional mapping of values to an index, the other a function call (that has no knowledge of any indexes on tables). You could observe similar behaviour from ST_Buffer or any other function that alters a geometry in a WHERE clause containing a statement that would normally use an index automatically.

If you rewrite your query so that ST_PointOnSurface is first calculated in a CTE (which involves no extra processing time, as you are calculating it anyway in your WHERE clause), then the index will be seen by the planner.

WITH pointonsurface (geom) AS 
  (SELECT 
      ST_PointOnSurface(geom) as geom 
   FROM tbl_vast_importance)
SELECT 
     id, geom 
FROM tbl_vast_importance a, pointonsurface b 
WHERE ST_Intersects(a.geom, b.geom); 

In my experience, rewriting queries in this way, so that any functions are calculated outside or any WHERE clauses or sub queries, can lead to order of magnitude increases in query time (due to avoiding table scans).

Answer 2

You can force it by manually adding a && (bounding box intersection) condition in your WHERE clause:

SELECT  
  id,
  geom
FROM tbl_vast_importance a
  JOIN tbl_vast_importance b
    ON a.geom && b.geom AND ST_Intersects(ST_PointOnSurface(a.geom),b.geom)
WHERE a.id != b.id