Software Development

Extending PostgreSQL: Complex Number Data Type

A few years ago I discussed creating user-defined types using PL/java. (Introduction to PostgreSQL PL/Java, part 4: User Defined Types) Today I follow up on the comment that we would use a standard C-language PostgreSQL extension in practice.

What is an extension?

A PostgreSQL extension is nothing more than bundled SQL statements and an optional native library. SQL-only extensions contain nothing but SQL statements.

Extensions provide two major benefits over the standalone scripts. First, the statements are loaded and unloaded as a unit. It is not possible to alter anything defined in an extension. This is a big win for maintainability and security. Second, extensions are versioned. This allows the SQL artifacts to be cleanly updated over time, something particularly important with user-defined types since simply dropping the UDT will cause the loss of data.

PostgreSQL has well-defined support for extensions. See Packaging Related Objects into an Extension.

Loading an extension

CREATE EXTENSION IF NOT EXISTS pg_complex;

Unloading an extension

DROP EXTENSION pg_complex;

It is important to remember that many if not most hosted database-as-a-service (DAAS) providers, e.g., Amazon RDS, will not allow you to install arbitrary extensions. You can still load SQL-only extensions by running the creation scripts manually but C-language extensions could cause problems in the future.

What are user-defined types (UDT)?

User-defined types are extensions to the database to support new types of data. SQL purists want tables to contain nothing but the standard primitives. Data that has a structure can be captured in a separate table and a foreign key, e.g., pulling an Address table out of a Customer table.

Many developers recognize that some data is more tightly bound than others and useful operations typically require more than one primitive value. Geographic position (latitude, longitude) and cryptographic information are classic examples. These objects are often so small that it does not make sense to have a separate table for them – compare a mailing address vs. the lat-long of that address.

(There may be other reasons to store the UDTs in a separate table, e.g., for improved security of sensitive information.)

Another use for UDTs is to add type safety to BLOBs. For instance it is reasonable to have user-defined functions that return the height or width of an image or the number of pages of a PDF document. You can easily write functions that accept BLOBs (bytea) but you can’t ensure that the value passed to the function is the appropriate type. Defining a UDT, e.g. pdf or jpeg, gives the developer a powerful tool.

In conclusion UDTs should be considered when 1) the object is meaningless if any element is missing or 2) the object would otherwise be a BLOB and you want to provide type-safe stored procedures and user-defined functions. Otherwise we should stick with the standard primitives.

Complex Numbers?

I will be creating a UDT for complex numbers below. I do not see a great unmet need to store complex numbers in relational databases but it is a good choice for educational purposes since it requires custom types, functions, casts, operators, and aggregate functions. The only thing missing is ordering.

This implementation uses a PostgreSQL composite type with a combination of SQL stored procedures and C user-defined functions. Complex numbers will be shown as (a, b) instead of the conventional a + bi but the latter is possible with additional work.

SQL definitions

This is a PostgreSQL extension so we should start by defining what we expect to see in our improved database.

Defining the complex UDT

The complex UDT consists of two fields – a real (re) component and an imaginary (im) component.

CREATE TYPE complex AS (re float8, im float8);

A very simple demonstration of its use follows.

$> CREATE TABLE t (c complex);
CREATE TABLE

-- insert a value. Note that we are inserting '(1,2)', not '1,2'.
$> INSERT INTO t VALUES((1,2));
INSERT 0 1

-- select full UDT
$> SELECT c FROM t;
   c   
-------
 (1,2)

-- select components. Note that we must surround field with parentheses.
$> SELECT (c).re, (c).im FROM t;
 re | im 
----+----
  1 |  2

Autopromoting floats to a complex numbers

It is easy to extract the real component of a complex number but still a pain to convert a real number to a complex number. PostgreSQL can do this transparently if we define a CAST.

CREATE OR REPLACE FUNCTION pgx_complex_from_int(int) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_from_bigint(bigint) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_from_numeric(numeric) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_from_numeric(float8) RETURNS complex AS $$
   SELECT ROW($1, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

Defining arithmetic operators

Complex numbers are numbers so we want to implement the standard arithmetic operators for them. All of the C function names are smurfed since they must unique. The SQL function names do not have to be unique since the function signature is also considered.

CREATE OPERATOR = (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_eq,
   NEGATOR = <>,
   HASHES,
   MERGES
);

CREATE OPERATOR  (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_ne,
   NEGATOR = <>,
   HASHES,
   MERGES
);

CREATE OPERATOR ~= (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_near,
   NEGATOR = <~>
);

CREATE OPERATOR  (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_not_near,
   NEGATOR = <>
);

CREATE OPERATOR - (
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_negate,
   NEGATOR = -
);

CREATE OPERATOR ~ (
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_conjugate,
   NEGATOR = ~
);

-- variants mixing 'complex' and 'numeric' types elided
CREATE OPERATOR + (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_add
);

-- variants mixing 'complex' and 'numeric' types elided
CREATE OPERATOR - (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgxdefine_complex_subtract
);

-- variants mixing 'complex' and 'numeric' types elided
CREATE OPERATOR * (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_multiply
);

-- variants mixing 'complex' and 'numeric' types elided
CREATE OPERATOR / (
   LEFT_ARG = complex,
   RIGHT_ARG = complex,
   PROCEDURE = pgx_complex_divide
);

Defining aggregate functions

Aggregate functions are one of the main ways to put intelligence into the database instead of treating it like a glorified file system. These are functions that operate on a collection of values and return an aggregate value of some type. Aggregate functions can have multiple parameters.

Don’t assume aggregates can only consume and produce numeric data. Consider a function that takes (x, y) pairs and produces a .png plot of the results.

Aggregate functions can optionally support window functions (http://www.postgresql.org/docs/9.4/static/tutorial-window.html). These are a relatively new feature and incredibly powerful. Implementation is conceptually simple – we only need functions to add or remove a value from the aggregate function’s ‘state’ – but in practice operations may be irreversible. That is the case here – if we compute ‘1e20 + 1 – 1e20′ the results should be ‘1’ but may be zero due to limited resolution.

--
-- UDT that keeps track of sums and sums-of-squares of a collection
-- of complex values. Tracking all three values allows us to compute
-- a wide variety of statistical values.
--
CREATE TYPE complex_accum AS (
   cnt  int,
   sum  complex,
   sofs complex
);

--
-- Calculate sum of a collection of complex values
--
CREATE AGGREGATE sum(complex) (
   sfunc = pg_complex_add,
   stype = complex_accum,
   initcond = '(0, 0)',
   -- msfunc = pg_complex_add,
   -- minvfunc = pg_complex_subtract,
   -- mstype = complex_accum,
   -- minitcond = (0, 0)'
);

--
-- Calculate average of a collection of complex values.
--
CREATE AGGREGATE avg(complex) (
   sfunc = pg_complex_accum,
   stype = complex_accum,
   finalfunc = pg_complex_avg
   -- msfunc = pg_complex_accum,
   -- minvfunc = pg_complex_disaccum,
   -- mstype = complex_accum,
);

(See: http://www.postgresql.org/docs/9.4/static/xaggr.html.)

Defining user-defined functions

We now know the functions and signatures that we must implement. In this case we can do most of the functions in pure SQL but choose to do a few in C to demonstrate advanced techniques.

Note: under TDD principles we should only implement enough to allow the tests to run. In this case the functions should return null. I’m not doing that here since the functions are so simple they can be verified at a glance. Any multi-line function should follow TDD principles and return null.

--
-- create functions implemented in C.
--
CREATE OR REPLACE FUNCTION pgx_complex_near(complex, complex)
RETURNS bool
AS 'pg_complex', 'pgx_complex_near'
LANGUAGE C IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_divide(complex, complex)
RETURNS complex
AS 'pg_complex', 'pgx_complex_divide'
LANGUAGE C IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION norm(complex)
RETURNS complex
AS 'pg_complex', 'pgx_complex_norm'
LANGUAGE C IMMUTABLE STRICT;

--
-- create functions implemented in SQL.
--
CREATE OR REPLACE FUNCTION pgx_complex_from_int(int) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_from_bigint(bigint) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_from_numeric(numeric) RETURNS complex AS $$
   SELECT ROW($1::float8, 0)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_eq(complex, complex) RETURNS bool AS $$
   SELECT $1.re = $2.re AND $1.im = $2.im;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_ne(complex, complex) RETURNS bool AS $$
   SELECT $1.re <> $2.re OR $1.im <> $2.im;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_not_near(complex, complex) RETURNS bool AS $$
   SELECT NOT pgx_complex_near($1, $2);
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_negate(complex) RETURNS complex AS $$
   SELECT ROW(-$1.re, -$1.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_conjugate(complex) RETURNS complex AS $$
   SELECT ROW($1.re, -$1.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_add(complex, complex) RETURNS complex AS $$
   SELECT ROW($1.re + $2.re, $1.im + $2.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_add_f8(float8, complex) RETURNS complex AS $$
   SELECT ROW($1 + $2.re, $2.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_add_f8(complex, float8) RETURNS complex AS $$
   SELECT ROW($1.re + $2, $1.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_subtract(complex, complex) RETURNS complex AS $$
   SELECT ROW($1.re - $2.re, $1.im - $2.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_subtract_f8(float8, complex) RETURNS complex AS $$
   SELECT ROW($1 - $2.re, -$2.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_subtract_f8(complex, float8) RETURNS complex AS $$
   SELECT ROW($1.re - $2, $1.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_multiply(complex, complex) RETURNS complex AS $$
   SELECT ROW($1.re * $2.re - $1.im * $2.im, $1.re * $2.im + $1.im * $2.re)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_multiply_f8(float8, complex) RETURNS complex AS $$
   SELECT ROW($1 * $2.re, $1 * $2.im)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION pgx_complex_multiply_f8(complex, float8) RETURNS complex AS $$
   SELECT ROW($1.re * $2, $1.im * $2)::complex;
$$ LANGUAGE SQL IMMUTABLE STRICT;

CREATE OR REPLACE FUNCTION magnitude(complex) RETURNS float8 AS $$
   SELECT sqrt($1.re * $1.re + $1.im * $1.im);
$$ LANGUAGE SQL IMMUTABLE STRICT;

Creating a skeleton extension

We’re now ready to create a skeleton extension. This extension follows test-driven development (TDD) practices – we want the minimal amount of code that fails. In this case that means the extension will load and define the user-defined functions and type but all functions and operators immediately return NULL.

The easiest way to do this is the PGXN utilities.

First, make sure that the following packages are installed:

  • pgxnclient
  • postgresql-server-dev-9.4
  • make
  • ruby
  • ruby2.1-dev
  • gcc

(This is for PostgreSQL 9.4 under Ubuntu. Adjust accordingly.)

Second, clone the github repository guedes/pgxn-utils.

Third, install these tools.

$ sudo pgxnclient install pgxn_utils

# verify utilities have been installed.
$ pgxn-utils help
PGXN Utils version: 0.1.4
Commands:
  pgxn-utils bundle [extension_name]  # Bundles the extension in a zip file
  pgxn-utils change [extension_name]  # Changes META's attributes in current extension
  pgxn-utils help [COMMAND]           # Describe available commands or one specific command
  pgxn-utils release filename         # Release an extension to PGXN
  pgxn-utils skeleton extension_name  # Creates an extension skeleton in current directory

Fourth, create a skeleton for a C-based PostgreSQL extension using our new utilities.

$ pgxn skeleton -m "Bear Giles <bgiles@coyotesong.com>" --template=c pg_complex
      create  pg_complex
      create  pg_complex/pg_complex.control
      create  pg_complex/.gitignore
      create  pg_complex/.template
      create  pg_complex/META.json
      create  pg_complex/Makefile
      create  pg_complex/README.md
      create  pg_complex/doc/pg_complex.md
      create  pg_complex/sql/pg_complex.sql
      create  pg_complex/sql/uninstall_pg_complex.sql
      create  pg_complex/src/pg_complex.c
      create  pg_complex/test/expected/base.out
      create  pg_complex/test/sql/base.sql

Fifth, edit the META.json, README.md and doc/pg_complex.md files to describe the extension. This would also be a good time to copy a LICENSE file into this directory if you have any plans to release the extension to others. Your future self will thank you for this documentation.

The META.json file allows us to specify dependencies among extensions.

Sixth, create a dummy implementation of each function that immediately returns null.

    #include "postgres.h"
    #include "fmgr.h"
     
    PG_MODULE_MAGIC;
     
    /**
     * Are two points "near" each other. This function requires reading
     * composite types.
     */
    PG_FUNCTION_INFO_V1(pgx_complex_near);
     
    Datum
    pgx_complex_near(PG_FUNCTION_ARGS) {
        PG_RETURN_NULL();
    }
     
    /**
     * Divide one complex number by another. This function requires reading
     * and returning composite types.
     */
    PG_FUNCTION_INFO_V1(pgx_complex_divide);
     
    Datum
    pgx_complex_divide(PG_FUNCTION_ARGS) {
        PG_RETURN_NULL();
    }
     
    /**
     * Scale a complex number so it on the unit circle. This function requires
     * reading and returning composite types.
     */
    PG_FUNCTION_INFO_V1(pgx_complex_norm);
     
    Datum
    pgx_complex_norm(PG_FUNCTION_ARGS) {
        PG_RETURN_NULL();
    }

Our problem is simple enough that the SQL stored procedures implemented the required functionality. More complex stored procedures could be implemented as plpsql stored procedures that return null.

Seventh, build the system.

$ make
$ sudo make install

You may need to load the extension before you call ‘make install’ the first time. It is not necessary to reload it afterwards.

$ sudo pgxn load ./

Eigth, run the tests.

The standard skeleton has support for regression tests.

$ make installcheck

The regression tests run all scripts under test/sql and verifies the results match the corresponding files in test/expected. The actual results are saved in results so it’s easy to write tests, modify the code as required, and then copy the file from results to test/expected once the desired behavior is seen.

You can also run tests through pgxn.

$ pgxn check -d somedb pg_complex

The optional pgTAP (http://pgtap.org/) extension gives us the ability to write xJunit-type tests. These tests will probably be more comfortable for developers than the regression tests.

For information integrating pgTAP into the build process see http://pgtap.org/integration.html and https://gkoenig.wordpress.com/2011/03/04/pgtap-unit-tests-for-postgresql/.

Ninth, install and deploy the extension outside of the test framework.

$ pgxn install --sudo -- pg_complex
$ pgxn load -d somedb --sudo -- pg_complex

You can undeploy and uninstall the extension using the analogous commands.

$ pgxn unload --sudo -- pg_complex
$ pgxn uninstall --sudo -- pg_complex

Tenth, publish the extension. You probably don’t want to do this with the skeleton implementation but this is a natural place to document the process. If we are a member of PGXN and wish to make our extension public we start by bundling our extension

$ pgxn bundle

and then upload it to https://manager.pgxn.org/.

Testing

As good test-driven development developers we start by writing our tests. In our case it’s straight SQL, or more precisely files that can be run through psql.

A typical test script is:

test/sql/math.sql

\set ECHO None
BEGIN;
\i sql/complex.sql
\set ECHO all

\set c1 (1,2)::complex
\set c2 (1,1)::complex
\set c3 (3,4)::complex
\set c4 (3,8)::complex

SELECT 1::complex AS a, (1::int8)::complex AS b, 1.0::complex AS c;

SELECT (1,2)::complex AS a, -(1,2)::complex AS b, ~(1,2)::complex AS c;

SELECT :c1 + (3,4)::complex AS a, 3 + :c1 AS b, :c1 + 3 AS c;

SELECT :c1 - (3,6)::complex AS a, 3 - :c1 AS b, :c1 - 3 AS c;

SELECT :c1 * (3,5)::complex AS a, 3 * :c1 AS b, :c1 * 3 AS c;
SELECT :c1 * (3,5)::complex AS a, 3.0::double precision * :c1 AS b, :c1 * 3.0 AS c;

SELECT :c4 / :c1  AS a, (:c4 / :c1) * :c1 = :c4 AS b;
SELECT :c4 / (2,0)::complex AS a, (2,0)::complex * (:c4 / (2,0)::complex)  = :c4 AS b;
SELECT :c4 / (0,2)::complex AS a, (0,2)::complex * (:c4 / (0,2)::complex) = :c4 AS b;
SELECT :c4 / 3 AS a, 3 * (:c4 / 3) = :c4 AS b;
SELECT 3 / :c4 AS a, :c4 * (3 / :c4) = 3::complex AS b;

--
-- check magnitude
--
SELECT magnitude(:c1) AS magnitude;
SELECT magnitude(:c2) AS magnitude;
SELECT magnitude(:c3) AS magnitude;

ROLLBACK;

The corresponding expected results are:

test/expected/math.out

\set ECHO None
\set c1 (1,2)::complex
\set c2 (1,1)::complex
\set c3 (3,4)::complex
\set c4 (3,8)::complex
SELECT 1::complex AS a, (1::int8)::complex AS b, 1.0::complex AS c;
   a   |   b   |   c   
-------+-------+-------
 (1,0) | (1,0) | (1,0)
(1 row)

SELECT (1,2)::complex AS a, -(1,2)::complex AS b, ~(1,2)::complex AS c;
   a   |    b    |   c    
-------+---------+--------
 (1,2) | (-1,-2) | (1,-2)
(1 row)

SELECT :c1 + (3,4)::complex AS a, 3 + :c1 AS b, :c1 + 3 AS c;
   a   |   b   |   c   
-------+-------+-------
 (4,6) | (4,2) | (4,2)
(1 row)

SELECT :c1 - (3,6)::complex AS a, 3 - :c1 AS b, :c1 - 3 AS c;
    a    |   b    |   c    
---------+--------+--------
 (-2,-4) | (2,-2) | (-2,2)
(1 row)

SELECT :c1 * (3,5)::complex AS a, 3 * :c1 AS b, :c1 * 3 AS c;
    a    |   b   |   c   
---------+-------+-------
 (-7,11) | (3,6) | (3,6)
(1 row)

SELECT :c1 * (3,5)::complex AS a, 3.0::double precision * :c1 AS b, :c1 * 3.0 AS c;
    a    |   b   |   c   
---------+-------+-------
 (-7,11) | (3,6) | (3,6)
(1 row)

SELECT :c4 / :c1  AS a, (:c4 / :c1) * :c1 = :c4 AS b;
     a     | b 
-----------+---
 (3.8,0.4) | t
(1 row)

SELECT :c4 / (2,0)::complex AS a, (2,0)::complex * (:c4 / (2,0)::complex)  = :c4 AS b;
    a    | b 
---------+---
 (1.5,4) | t
(1 row)

SELECT :c4 / (0,2)::complex AS a, (0,2)::complex * (:c4 / (0,2)::complex) = :c4 AS b;
    a     | b 
----------+---
 (4,-1.5) | t
(1 row)

SELECT :c4 / 3 AS a, 3 * (:c4 / 3) = :c4 AS b;
          a           | b 
----------------------+---
 (1,2.66666666666667) | t
(1 row)

SELECT 3 / :c4 AS a, :c4 * (3 / :c4) = 3::complex AS b;
                   a                    | b 
----------------------------------------+---
 (0.123287671232877,-0.328767123287671) | t
(1 row)

--
-- check magnitude
--
SELECT magnitude(:c1) AS magnitude;
    magnitude     
------------------
 2.23606797749979
(1 row)

SELECT magnitude(:c2) AS magnitude;
    magnitude    
-----------------
 1.4142135623731
(1 row)

SELECT magnitude(:c3) AS magnitude;
 magnitude 
-----------
         5
(1 row)

ROLLBACK;

It is probably easiest to create the ‘expected’ file by running the test once, getting the results from results/math.out, and editing that file to show the expected results. In a pure TDD implementation all of the tests should initially return null but we’ve already defined many of the functions above.

Implementation

There are three SQL-language functions to add. Details matter – the sum of no values is well-defined as 0 + 0i but the average of no values is undefined (null), not any particular value.

--
-- accumulator function is similar to float8_accum. Question: should the result
-- be the product of p * p or the product of p * ~p ?
--
CREATE OR REPLACE FUNCTION pgx_complex_accum(complex_accum, complex) RETURNS complex_accum AS $$
   SELECT CASE WHEN $1 IS NULL THEN 1 ELSE $1.cnt + 1 END,
          CASE WHEN $1 IS NULL THEN $2 ELSE $1.sum + $2 END,
          CASE WHEN $1 IS NULL THEN $2 * ~$2 ELSE $1.sofs + $2 * ~$2 END;
$$ LANGUAGE SQL;

--
-- disaccumulator(?) function is similar to pgx_complex_accum. It is required in order
-- to implement windowing functions.
--
CREATE OR REPLACE FUNCTION pgx_complex_disaccum(complex_accum, complex) RETURNS complex_accum AS $$
   SELECT pgx_complex_accum($1, -$2);
$$ LANGUAGE SQL;

--
-- average function returns quotient of sum over count.
--
CREATE OR REPLACE FUNCTION pgx_complex_avg(complex_accum) RETURNS complex AS $$
   SELECT CASE WHEN $1 IS NULL THEN NULL
               WHEN $1.cnt = 0 THEN (0,0)::complex
               ELSE $1.sum / $1.cnt END;
$$ LANGUAGE SQL;

The first C-language function demonstrates how to read a composite value and return a primitive. In this case we get the ‘re’ and ‘im’ components by name.

PG_MODULE_MAGIC;

/**
 * Test complex numbers for proximity. This avoids the problems with testing floats
 * and doubles but does not guarantee absolute equality.
 */
PG_FUNCTION_INFO_V1(pgx_complex_near);

Datum
pgx_complex_near(PG_FUNCTION_ARGS) {
    double re[2];
    double im[2];
    double p, q;
    int i;

    // unwrap values.    
    for (i = 0; i < 2; i++) {
        HeapTupleHeader t = PG_GETARG_HEAPTUPLEHEADER(i);
        bool isnull[2];

        Datum dr = GetAttributeByName(t, "re", &isnull[0]);
        Datum di = GetAttributeByName(t, "im", &isnull[1]);

        // STRICT prevents the 'complex' value from being null but does
        // not prevent its components from being null.        
        if (isnull[0] || isnull[1]) {
            PG_RETURN_NULL();
        }
        
        re[i] = DatumGetFloat8(dr);
        im[i] = DatumGetFloat8(di);
    }

    // compute distance between points, distance of points from origin.
    p = hypot(re[0] - re[1], im[0] - im[1]);
    q = hypot(re[0], im[0]) + hypot(re[1], im[1]);
    
    if (q == 0) {
        PG_RETURN_BOOL(1);
    }
    
    // we consider the points 'near' each other if the distance between them is small
    // relative to the size of them. 
    PG_RETURN_BOOL(p / q < 1e-8); 
}

The second case returns a composite value. There are two ways to return composite values. This is the older way and requires a little more work. The newer way requires everything be returned as a string – this has a modest cost with primitive values but it could be costly to marshal and unmarshal user-defined types.

/**
 * Divide complex number by another. We do this by multiplying nominator and denominator
 * by the conjugate of the denominator. The denominator then becomes the scalar square of
 * the magnitude of the number.
 */
PG_FUNCTION_INFO_V1(pgx_complex_divide);

Datum
pgx_complex_divide(PG_FUNCTION_ARGS) {
    TupleDesc tupdesc;
    HeapTuple tuple;
    double re[2];
    double im[2];
    int i;
    double q;
    Datum datum[2];
    bool isnull[2];
 
    // build a tuple descriptor for our result type 
    if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE)
        ereport(ERROR,
                (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
                 errmsg("function returning record called in context "
                        "that cannot accept type record")));

    // unwrap values.    
    for (i = 0; i < 2; i++) {
        HeapTupleHeader t = PG_GETARG_HEAPTUPLEHEADER(i);
        bool isnull[2];
        Datum dr, di;

        dr = GetAttributeByName(t, "re", &isnull[0]);
        di = GetAttributeByName(t, "im", &isnull[1]);

        // STRICT prevents the 'complex' value from being null but does
        // not prevent its components from being null.        
        if (isnull[0] || isnull[1]) {
            PG_RETURN_NULL();
        }
        
        re[i] = DatumGetFloat8(dr);
        im[i] = DatumGetFloat8(di);
    }

    // the denominator is the square of the magnitude of the divisor.
    q = re[1] * re[1] + im[1] * im[1];
    
    // should I throw error instead of returning null?
    if (q == 0.0) {
        PG_RETURN_NULL();
    }

    datum[0] = Float8GetDatum((re[0] * re[1] + im[0] * im[1]) / q);
    datum[1] = Float8GetDatum((im[0] * re[1] - im[1] * re[0]) / q);

    BlessTupleDesc(tupdesc);
    tuple = heap_form_tuple(tupdesc, datum, isnull);
 
    PG_RETURN_DATUM(HeapTupleGetDatum(tuple));
}

The third example also consumes and produces a composite value.

/**
 * Calculate the norm of a complex number. This is the complex number on the unit
 * circle so that magnitude(norm(x)) = 1 and magnitude(x) * norm(x) = x.
 */
PG_FUNCTION_INFO_V1(pgx_complex_norm);

Datum
pgx_complex_norm(PG_FUNCTION_ARGS) {
    HeapTupleHeader t = PG_GETARG_HEAPTUPLEHEADER(0);
    TupleDesc tupdesc;
    HeapTuple tuple;
    double re;
    double im;
    bool isnull[2];
    Datum datum[2];
    double m;
 
    // build a tuple descriptor for our result type 
    if (get_call_result_type(fcinfo, NULL, &tupdesc) != TYPEFUNC_COMPOSITE)
        ereport(ERROR,
                (errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
                 errmsg("function returning record called in context "
                        "that cannot accept type record")));
        
    // unwrap values.    
    datum[0] = GetAttributeByName(t, "re", &isnull[0]);
    datum[1] = GetAttributeByName(t, "im", &isnull[1]);

    // STRICT prevents the 'complex' value from being null but does
    // not prevent its components from being null.        
    if (isnull[0] || isnull[1]) {
        PG_RETURN_NULL();
    }
        
    re = DatumGetFloat8(datum[0]);
    im = DatumGetFloat8(datum[1]);

    m = hypot(re, im);
   
    // should I throw error instead of returning null?
    if (m == 0.0) {
        PG_RETURN_NULL();
    } 

    datum[0] = Float8GetDatum(re / m);
    datum[1] = Float8GetDatum(im / m);

    BlessTupleDesc(tupdesc);
    tuple = heap_form_tuple(tupdesc, datum, isnull);
 
    PG_RETURN_DATUM(HeapTupleGetDatum(tuple));
}

Wrap up

This wraps up the complex number extension. Most extensions will be far more simple but this problem was chosen precisely because it demonstrates how deeply you can integrate an extension.

Source code

  • http://github.com/beargiles/pg-complex
  • http://pgxn.org/dist/complex
  • Additional Resources

Subscribe
Notify of
guest

This site uses Akismet to reduce spam. Learn how your comment data is processed.

0 Comments
Inline Feedbacks
View all comments
Back to top button