NRPP, Averter of Disaster

Far and away the most amusing part of the NRPP Business Case that the government quietly placed online this March is the “Financial Benefits & Lifecycle Costs” section, the part where they are supposed to make the dollars-and-cents case for the project.

In order to ensure transparency and instill public confidence in the project, all the dollar figures have been replaced with “$X” in the text. So awesome!

But just as good, in the Lower Risks Exposure section, they outline in clear prose some of the resource management disasters of the past 10 years:

  • A dam collapsed and caused significant property damage after the NRS failed to enforce the need for repairs and maintenance. Contributing factors included a lack of effective record keeping, inaccurate risk/consequence ratings given to the dam, and a lack of action to address decades of warnings that the dam was in poor condition. Settlements paid by government to affected parties totaled more than $X M
  • A District Manager approved a forest development plan and awarded multiple cut blocks in a watershed; however, the decision was inappropriate and not durable because the plan was contrary to regulations and it did not consider impacts to endangered animal populations. The case was settled for $X M
  • An oil and gas company acquired a tenure in northeastern British Columbia only to find out that the NRS had not disclosed that the tenure was in an area that was considered historically and spiritually significant by a local aboriginal group. The matter proceeded through an extensive litigation process and settled for $X M

That’s some pretty bad stuff! Not the kind of thing that makes it into the average government press release!

But, being information technology people, they cannot help themselves, they look at this bad stuff and they think, “if only people had the right information at the right time, none of this would have happened”.

NRPP, Averter of Disaster

Through the implementation of an integrated spatial and operational database, a risk-based framework for compliance and enforcement activities across the province, and through facilitating information sharing across lines of business, NRPP will contribute to a significant reduction in legal risk and potentially a reduction in litigation and settlement costs.

Aren’t IT people cute! Don’t you just want to stroke their fur and take them home?!?

Non-IT people might have a different reaction to that list of calamities. They might think, “political pressures that make speed of permit approval the primary metric of good resource management are causing corners to be cut and very expensive/damaging mistakes to be made.”

Unfortunately, those are the kinds of thoughts that will quickly end careers if expressed out loud in government.

So instead, we have a new $57M IT megaproject. The kind of project that has the maximum likelihood of failure and delay: a multi-year, multi-million-dollar, maximum scope project.

Don’t worry though, it will be delivered “on time and on budget”.

NRPP, On Time and On Budget

As Vaughn Palmer said last week,

“Another day, another bogus claim from the B.C. Liberals that a major capital project is on time and on budget.”

We in the IT consulting field are perhaps numb to the eternal cycle of cutting scope, pushing the delivery date, and declaring victory, but in the cut-throat world of politics, these things have consequences.

At the very least, you’ll be publicly mocked when you do it, as on Tuesday the NDP mocked the Liberals for repeated non-delivery of an online Civil Resolution Tribunal.

NRPP, On Time and On Budget

Mind you, most of these wounds are self-inflicted, as this government just will not let go of the core principles of enterprise IT failure:

  • make the budget as large as possible;
  • make the scope as broad as possible; and,
  • hire international IT consultancies so large that they don’t care whether your project succeeds or not, just so long as they can bill out as much as possible.

These are the guiding principles that have brought them ICM, MyEdBC, Panorama, and enough other miscues that the Vancouver Sun felt compelled to run a whole series of articles on IT failure.

And there’s another one in progress!

The $57.2 million Natural Resource Permitting Project (NRPP) Phase One is now at least a year behind schedule, based on a business plan document that was quietly placed online this March.

This is the same NRPP business plan I requested via FOI in 2014, and was given 200 blank pages of PDF in reply. The plan released March 2016 still has all the dollar figures blanked out. That’s right, they released a “business plan” without the numbers!

According to the “Investment Roadmap” year one of the project, fiscal 2014/15 which ended 12 months ago, should have delivered the following results:

  • New legislative framework in place
  • High-priority business processes re-engineered
  • Initial foundational integrated spatial-temporal operational database available
  • IM/IT and priority finance operating model and processes re-engineered
  • New Common Client architecture available
  • Transformation governance model in place

Has any of this happened?

The project capital plan, written last fall, lists similar capabilities as items to be delivered in the future.

According to the business plan, on the basis of which $57M capital dollars were committed in 2014, the project is now at least a year behind.

As we can see from the discrepency between the older business plan and the capital plan, the usual process has begun: the goal posts are on the move. Next up, the bogus claims.

Fortunately, the international consultants are billing hard so we’ll be back on schedule shortly. In the end everything will be delivered “on time, and on budget”, I guarantee it.

OGR FDW Update

I’ve had a productive couple of weeks here, despite the intermittently lovely weather and the beginning of Little League baseball season (not coaching, just supporting my pitcher-in-training).

13 Days

The focus of my energies has been a long-awaited (by me) update to the OGR FDW extension for PostgreSQL. By binding the multi-format OGR library into PostgreSQL, we get access to the many formats supported by OGR, all with just one piece of extension code.

As usual, the hardest part of the coding was remembering how things worked in the first place! But after getting my head back in the game the new code flowed out and now I can reveal the new improved OGR FDW!

OGR FDW Update

The new features are:

  • Column name mapping between OGR layers and PgSQL tables is now completely configurable. The extension will attempt to guess mapping automagically, using names and type consistency, but you can over-ride mappings using the table-level column_name option.
  • Foreign tables are now updateable! That means, for OGR sources that support it, you can run INSERT, UPDATE and DELETE commands on your OGR FDW tables and the changes will be applied to the source.

    • You can control which tables and foreign servers are updateable by setting the UPDATEABLE option on the foreign server and foreign table definitions.
  • PostgreSQL 9.6 is supported. It’s not released yet, but we can now build against it.
  • Geometry type and spatial reference system are propogated from OGR. If your OGR source defines a geometry type and spatial reference identifier, the FDW tables in PostGIS will now reflect that, for easier integration with your local geometry data.
  • GDAL2 and GDAL1 are supported. Use of GDAL2 syntax has been made the default in the code-base, with mappings back to GDAL1 for compatibility, so the code is now future-ready.
  • Regression tests and continuous integration are in place, for improved code reliability. Thanks to help from Even Roualt, we are now using Travis-CI for integration testing, and I’ve enabled a growing number of integration tests.

As usual, I’m in debt to Regina Obe for her timely feedback and willingness to torture-test very fresh code.

For now, early adopters can get the code by cloning and building the project master branch, but I will be releasing a numbered version in a week or two when any obvious bugs have been shaken out.

Parallel PostGIS Joins

In my earlier post on new parallel query support coming in PostgreSQL 9.6 I was unable to come up with a parallel join query, despite much kicking and thumping the configuration and query.

Parallel PostGIS Joins

It turns out, I didn’t have all the components of my query marked as PARALLEL SAFE, which is required for the planner to attempt a parallel plan. My query was this:

 SELECT Count(*) 
  FROM pd 
  JOIN pts 
  ON pd.geom && pts.geom
  AND _ST_Intersects(pd.geom, pts.geom);

And _ST_Intersects() was marked as safe, but I neglected to mark the function behind the && operator – geometry_overlaps – as safe. With both functions marked as safe, and assigned a hefty function cost of 1000, I get this query:

 Nested Loop  
 (cost=0.28..1264886.46 rows=21097041 width=2552) 
 (actual time=0.119..13876.668 rows=69534 loops=1)
   ->  Seq Scan on pd  
   (cost=0.00..14271.34 rows=69534 width=2512) 
   (actual time=0.018..89.653 rows=69534 loops=1)
   ->  Index Scan using pts_gix on pts  
   (cost=0.28..17.97 rows=2 width=40) 
   (actual time=0.147..0.190 rows=1 loops=69534)
         Index Cond: (pd.geom && geom)
         Filter: _st_intersects(pd.geom, geom)
         Rows Removed by Filter: 2
 Planning time: 8.365 ms
 Execution time: 13885.837 ms

Hey wait! That’s not parallel either!

It turns out that parallel query involves a secret configuration sauce, just like parallel sequence scan and parellel aggregate, and naturally it’s different from the other modes (gah!)

The default parallel_tuple_cost is 0.1. If we reduce that by an order of magnitude, to 0.01, we get this plan instead:

 (cost=1000.28..629194.94 rows=21097041 width=2552) 
 (actual time=0.950..6931.224 rows=69534 loops=1)
   Number of Workers: 3
   ->  Nested Loop  
   (cost=0.28..628194.94 rows=21097041 width=2552) 
   (actual time=0.303..6675.184 rows=17384 loops=4)
         ->  Parallel Seq Scan on pd  
         (cost=0.00..13800.30 rows=22430 width=2512) 
         (actual time=0.045..46.769 rows=17384 loops=4)
         ->  Index Scan using pts_gix on pts  
         (cost=0.28..17.97 rows=2 width=40) 
         (actual time=0.285..0.366 rows=1 loops=69534)
               Index Cond: (pd.geom && geom)
               Filter: _st_intersects(pd.geom, geom)
               Rows Removed by Filter: 2
 Planning time: 8.469 ms
 Execution time: 6945.400 ms

Ta da! A parallel plan, and executing almost twice as fast, just like the doctor ordered.


Mostly the parallel support in core “just works” as advertised. PostGIS does need to mark our functions as quite costly, but that’s reasonable since they actually are quite costly. What is not good is the need to tweak the configuration once the functions are properly costed:

  • Having to cut parallel_tuple_cost by a factor of 10 for the join case is not any good. No amount of COST increases seemed to have an effect, only changing the core parameter did.
  • Having to increase the cost of functions used in aggregates by a factor of 100 over cost of functions used in sequence filters is also not any good.

So, with a few more changes to PostGIS, we are quite close, but the planner for parallel cases needs to make more rational use of function costs before we have achieved parallel processing nirvana for PostGIS.

Parallel PostGIS

Parallel query support in PostgreSQL in the upcoming 9.6 release will be available for a number of query types: sequence scans, aggregates and joins. Because PostGIS tends to involve CPU-intensive calculations on geometries, support for parallel query has been at the top of our request list to the core team for a long time. Now that it is finally arriving the question is: does it really help?

Parallel PostGIS


  • With some adjustments to function COST both parallel sequence scan and parallel aggregation deliver very good parallel performance results.
  • The cost adjustments for sequence scan and aggregate scan are not consistent in magnitude.
  • Parallel join does not seem to work for PostGIS indexes yet, but perhaps there is some magic to learn from PostgreSQL core on that.


In order to run these tests yourself, you will need to check out and build:

For testing, I used the set of 69534 polling divisions defined by Elections Canada.

shp2pgsql -s 3347 -I -D -W latin1 PD_A.shp pd | psql parallel

It’s worth noting that this data set is, in terms of number of rows very very small as databases go. This will become important as we explore the behaviour of the parallel processing, because the assumptions of the PostgreSQL developers about what constitutes a “parallelizable load” might not match our assumptions in the GIS world.

With the data loaded, we can do some tests on parallel query. Note that there are some new configuration options for parallel behaviour that will be useful during testing:

  • max_parallel_degree sets the maximum degree of parallelism for an individual parallel operation. Default 0.
  • parallel_tuple_cost sets the planner’s estimate of the cost of transferring a tuple from a parallel worker process to another process. The default is 0.1.
  • parallel_setup_cost sets the planner’s estimate of the cost of launching parallel worker processes. The default is 1000.
  • force_parallel_mode allows the use of parallel queries for testing purposes even in cases where no performance benefit is expected. Default ‘off’.

Parallel Sequence Scan

Before we can test parallelism, we need to turn it on! The default max_parallel_degree is zero, so we need a non-zero value. For my tests, I’m using a 2-core laptop, so:

SET max_parallel_degree=2;

Now we are ready to run a query with a spatial filter. Using EXPLAIN ANALYZE suppressed the actual answer in favour of returning the query plan and the observed execution time:

  SELECT Count(*) 
    FROM pd 
    WHERE ST_Area(geom) > 10000;

And the answer we get back is:

 (cost=14676.95..14676.97 rows=1 width=8) 
 (actual time=757.489..757.489 rows=1 loops=1)
   ->  Seq Scan on pd  
   (cost=0.00..14619.01 rows=23178 width=0) 
   (actual time=0.160..747.161 rows=62158 loops=1)
         Filter: (st_area(geom) > 10000)
         Rows Removed by Filter: 7376
 Planning time: 0.137 ms
 Execution time: 757.553 ms

Two things we can learn here:

  • There is no parallelism going on here, the query plan is just a single-threaded one.
  • The single-threaded execution time is about 750ms.

Now we have a number of options to fix this problem:

  • We can force parallelism using force_parallel_mode, or
  • We can force parallelism by decreasing the parallel_setup_cost, or
  • We can adjust the cost of ST_Area() to try and get the planner to do the right thing automatically.

It turns out that the current definition of ST_Area() has a default COST setting, so it is considered to be no more or less expensive than something like addition or substraction. Since calculating area involves multiple floating point operations per polygon segment, that’s a stupid cost.

In general, all PostGIS functions are going to have to be reviewed and costed to work better with parallelism.

If we redefine ST_Area() with a big juicy cost, things might get better.

  AS '$libdir/postgis-2.3','area'
  COST 100;

Now the query plan for our filter is much improved:

Finalize Aggregate  
(cost=20482.97..20482.98 rows=1 width=8) 
(actual time=345.855..345.856 rows=1 loops=1)
->  Gather  
   (cost=20482.65..20482.96 rows=3 width=8) 
   (actual time=345.674..345.846 rows=4 loops=1)
     Number of Workers: 3
     ->  Partial Aggregate  
         (cost=19482.65..19482.66 rows=1 width=8) 
         (actual time=336.663..336.664 rows=1 loops=4)
           ->  Parallel Seq Scan on pd  
               (cost=0.00..19463.96 rows=7477 width=0) 
               (actual time=0.154..331.815 rows=15540 loops=4)
                 Filter: (st_area(geom) > 10000)
                 Rows Removed by Filter: 1844
Planning time: 0.145 ms
Execution time: 349.345 ms

Three important things to note:

  • We have a parallel query plan!
  • Some of the execution results output are wrong! They say that only 1844 rows were removed by the filter, but in fact 7376 were (as we can confirm by running the queries without the EXPLAIN ANALYZE). This is a known limitation, reporting on the results of only one parallel worker, which (should) maybe, hopefully be fixed before 9.6 comes out.
  • The execution time has been halved, just as we would hope for a 2-core machine!

Now for the disappointing part, try this:

  SELECT ST_Area(geom) 
    FROM pd;

Even though the work being carried out (run ST_Area() on 70K polygons) is exactly the same as in our first example, the planner does not parallelize it, because the work is not in the filter.

 Seq Scan on pd  
 (cost=0.00..31654.84 rows=69534 width=8) 
 (actual time=0.130..722.286 rows=69534 loops=1)
 Planning time: 0.078 ms
 Execution time: 727.344 ms

For geospatial folks, who tend to do a fair amount of expensive calculation in the SELECT parameters, this is a bit disappointing. However, we still get impressive parallelism on the filter!

Parallel Aggregation

The aggregate most PostGIS users would like to see parallelized is ST_Union() so it’s worth explaining why that’s actually a little hard.

PostgreSQL Aggregates

All aggregate functions in PostgreSQL consist of at least two functions:

  • A “transfer function” that takes in a value and a transfer state, and adds the value to the state. For example, the Avg() aggregate has a transfer state consisting of the sum of all values seen so far, and the count of all values processed.
  • A “final function” that takes in a transfer state and converts it to the final aggregate output value. For example, the Avg() aggregate final function divides the sum of all values by the count of all values and returns that number.

For parallel processing, PostgreSQL adds a third kind of function:

  • A “combine” function, that takes in two transfer states and outputs a singe combined state. For the Avg() aggregate, this would add the sums from each state and counts from each state and return that as the new combined state.

So, in order to get parallel processing in an aggregate, we need to define “combine functions” for all the aggregates we want parallelized. That way the master process can take the completed transfer states of all parallel workers, combine them, and then hand that final state to the final function for output.

To sum up, in parallel aggregation:

  • Each worker runs “transfer functions” on the records it is responsible for, generating a partial “transfer state”.
  • The master takes all those partial “transfer states” and “combines” them into a “final state”.
  • The master then runs the “final function” on the “final state” to get the completed aggregate.

Note where the work occurs: the workers only run the transfer functions, and the master runs both the combine and final functions.

PostGIS ST_Union Aggregate

One of the things we are proud of in PostGIS is the performance of our ST_Union() implementation, which gains performance from the use of a cascaded union algorithm.

Cascaded union involves the following steps:

  • Collects all the geometries of interest into an array (aggregate transfer function), then
  • Builds a tree on those geometries and unions them from the leaves of the tree upwards (aggregate final function).

Note that all the hard work happens in the final step. The transfer functions (which is what would be run on the workers) do very little work, just gathering geometries into an array.

Converting this process into a parallel one by adding a combine function that does the union would not make things any faster, because the combine step also happens on the master. What we need is an approach that does more work during the transfer function step.

PostGIS ST_MemUnion Aggregate

“Fortunately” we have such an aggregate, the old union implementation from before we added “cascaded union”. The “memory friendly” union saves memory by not building up the array of geometries in memory, at the cost of spending lots of CPU unioning each input geometry into the transfer state.

In that respect, it is the perfect example to use for testing parallel aggregate.

The non-parallel definition of ST_MemUnion() is this:

  basetype = geometry,
  sfunc = ST_Union,
  stype = geometry

No special types or functions required: the transfer state is a geometry, and as each new value comes in the two-parameter version of the ST_Union() function is called to union it onto the state. There is no final function because the transfer state is the output value. Making the parallel version is as simple as adding a combine function that also uses ST_Union() to merge the partial states:

  basetype = geometry,
  sfunc = ST_Union,
  combinefunc = ST_Union,
  stype = geometry

Now we can run an aggregation using ST_MemUnion() to see the results. We will union the polling districts of just one riding, so 169 polygons:

  SELECT ST_Area(ST_MemUnion(geom)) 
    FROM pd 
    WHERE fed_num = 47005;

Hm, no parallelism in the plan, and an execution time of 3.7 seconds:

 (cost=14494.92..14495.18 rows=1 width=8) 
 (actual time=3784.781..3784.782 rows=1 loops=1)
   ->  Seq Scan on pd  
   (cost=0.00..14445.17 rows=199 width=2311) 
   (actual time=0.078..49.605 rows=169 loops=1)
         Filter: (fed_num = 47005)
         Rows Removed by Filter: 69365
 Planning time: 0.207 ms
 Execution time: 3784.997 ms

We have to bump the cost of the two parameter version of ST_Union() up to 10000 before parallelism kicks in:

CREATE OR REPLACE FUNCTION ST_Union(geom1 geometry, geom2 geometry)
  RETURNS geometry
  AS '$libdir/postgis-2.3','geomunion'
  COST 10000;

Now we get a parallel execution! And the time drops down substantially, though not quite a 50% reduction.

 Finalize Aggregate  
 (cost=16536.53..16536.79 rows=1 width=8) 
 (actual time=2263.638..2263.639 rows=1 loops=1)
   ->  Gather  
   (cost=16461.22..16461.53 rows=3 width=32) 
   (actual time=754.309..757.204 rows=4 loops=1)
         Number of Workers: 3
         ->  Partial Aggregate  
         (cost=15461.22..15461.23 rows=1 width=32) 
         (actual time=676.738..676.739 rows=1 loops=4)
               ->  Parallel Seq Scan on pd  
               (cost=0.00..13856.38 rows=64 width=2311) 
               (actual time=3.009..27.321 rows=42 loops=4)
                     Filter: (fed_num = 47005)
                     Rows Removed by Filter: 17341
 Planning time: 0.219 ms
 Execution time: 2264.684 ms

The punchline though, is what happens when we run the query using a single-threaded ST_Union() with cascaded union:

  SELECT ST_Area(ST_Union(geom)) 
    FROM pd 
    WHERE fed_num = 47005;

Good algorithms beat brute force still:

 (cost=14445.67..14445.93 rows=1 width=8) 
 (actual time=2031.230..2031.231 rows=1 loops=1)
   ->  Seq Scan on pd  
   (cost=0.00..14445.17 rows=199 width=2311) 
   (actual time=0.124..66.835 rows=169 loops=1)
         Filter: (fed_num = 47005)
         Rows Removed by Filter: 69365
 Planning time: 0.278 ms
 Execution time: 2031.887 ms

The open question is, “can we combine the subtlety of the cascaded union algorithm with the brute force of parallel execution”?

Maybe, but it seems to involve magic numbers: if the transfer function paused every N rows (magic number) and used cascaded union to combine the geometries received thus far, it could possibly milk performance from both smart evaluation and multiple CPUs. The use of a magic number is concerning however, and the approach would be very sensitive to the order in which rows arrived at the transfer functions.

Parallel Join

To test parallel join, we’ll build a synthetic set of points, such that each point falls into one polling division polygon:

  ST_PointOnSurface(geom)::Geometry(point, 3347) AS geom, 
  gid, fed_num 
FROM pd;

  ON pts USING GIST (geom);

Points and Polling Divisions

A simple join query looks like this:

 SELECT Count(*) 
  FROM pd 
  JOIN pts 
  ON ST_Intersects(pd.geom, pts.geom);

But the query plan has no parallel elements! Uh oh!

(cost=222468.56..222468.57 rows=1 width=8) 
(actual time=13830.361..13830.362 rows=1 loops=1)
   ->  Nested Loop  
   (cost=0.28..169725.95 rows=21097041 width=0) 
   (actual time=0.703..13815.008 rows=69534 loops=1)
         ->  Seq Scan on pd  
         (cost=0.00..14271.34 rows=69534 width=2311) 
         (actual time=0.086..90.498 rows=69534 loops=1)
         ->  Index Scan using pts_gix on pts  
         (cost=0.28..2.22 rows=2 width=32) 
         (actual time=0.146..0.189 rows=1 loops=69534)
               Index Cond: (pd.geom && geom)
               Filter: _st_intersects(pd.geom, geom)
               Rows Removed by Filter: 2
 Planning time: 6.348 ms
 Execution time: 13843.946 ms

The plan does involve a nested loop, so there should be an opportunity for parallel join to work magic. Unfortunately no variation of the query or the parallel configuration variables, or the function costs will change the situation: the query refuses to parallelize!

SET parallel_tuple_cost=0.001;
SET force_parallel_mode=on;
SET parallel_setup_cost=1;

The ST_Intersects() function is actually a SQL wrapper on top of the && operator and the _ST_Intersects() function, but unwrapping it and using the components directly also has no effect.

 SELECT Count(*) 
  FROM pd 
  JOIN pts 
  ON pd.geom && pts.geom
  AND _ST_Intersects(pd.geom, pts.geom);

The only variant I could get to parallelize omitted the && index operator.

 SELECT *        
  FROM pd 
  JOIN pts 
  ON _ST_Intersects(pd.geom, pts.geom);

Unfortunately without the index operator the query is so inefficient it doesn’t matter that it’s being run in parallel, it will take days to run to completion.

 (cost=1000.00..721919734.88 rows=1611658891 width=2552)
   Number of Workers: 2
   ->  Nested Loop  
   (cost=0.00..576869434.69 rows=1611658891 width=2552)
         Join Filter: _st_intersects(pd.geom, pts.geom)
         ->  Parallel Seq Scan on pd  
         (cost=0.00..13865.73 rows=28972 width=2512)
         ->  Seq Scan on pts  
         (cost=0.00..1275.34 rows=69534 width=40)

So, thus far, parallel query seems to be a wet squib for PostGIS, though I hope with some help from PostgreSQL core we can figure out where the problem lies.

UPDATE: Marking the geometry_overlaps function which is bound to the && operator as PARALLEL SAFE allows PostgreSQL to generate parallel join plans when the index is in place.


While it is tempting to think “yay, parallelism! all my queries will run $ncores times faster!” in fact parallelism still only applies in a limited number of cases:

  • When there is a sequence scan large (costly) enough to be worth parallelizing.
  • When there is an aggregate large (costly) enough to be worth parallelizing, and the aggregate function can actually parallize the work effectively.
  • (Theoretically) when there is a (nested loop) join large (costly) enough to be worth parallelizing.

Additionally there is still work to be done on PostGIS for optimal use of the parallel features we have available:

  • Every function is going to need a cost, and those costs may have to be quite high to signal to the planner that we are not garden variety computations.
  • Differences in COST adjustments for different modes need to be explored: why was a 10000 cost needed to kick the aggregation into action, while a 100 cost sufficed for sequence scan?
  • Aggregation functions that currently back-load work to the final function may have to be re-thought to do more work in the transfer stage.
  • Whatever issue is preventing our joins from parallelizing needs to be tracked down.

All that being said, the potential is to see a large number of queries get $ncores faster, so this promises to be the most important core development we’ve seen since the extension framework arrived back in PostgreSQL 9.1.