I'm trying to run this but after 10 hours I get the error "cannot have more than 2^32-1 subtransactions in a transaction". What is your advice for getting around this issue? I'm hesitant to split up this problem but I think if I split it at the 30° parallel maybe that will work... I'll try that since this is kind of an uncommon thing I think

I'm using postgresql 11

    CREATE TABLE cities_union AS SELECT k, COALESCE(ST_Difference(geom, (SELECT st_union(ST_MakeValid(ST_SnapToGrid(b.geom, 0.0001)))
    FROM test600 b
    WHERE safe_isect(a.geom, b.geom)
         and a.area > 0
         and b.area > 0
         AND a.area > b.area
         and wof_id is not null
         )
      ), 
      a.geom) 
    FROM test600 a 
    where a.area > 0 
    ;

I'm using this function:

    CREATE OR REPLACE FUNCTION safe_isect(geom_a geometry, geom_b geometry)
    RETURNS boolean AS
    $$
    BEGIN
        RETURN ST_Intersects(geom_a, geom_b);
        EXCEPTION
            WHEN OTHERS THEN
                BEGIN
                    RETURN ST_Intersects(ST_Buffer(geom_a, 0.0000001), ST_Buffer(geom_b, 0.0000001));
                    EXCEPTION
                        WHEN OTHERS THEN
                            RETURN ST_GeomFromText('POLYGON EMPTY');
        END;
    END
    $$
    LANGUAGE 'plpgsql' STABLE STRICT;

Previously when I ran the following command I was able to get to 30 hours of runtime but I had an error with one of the geometries which caused it to stop and is why I modified the ST_Difference part above

    CREATE TABLE cities_union AS SELECT k, COALESCE(ST_Difference(geom, (SELECT st_union(b.geom, 0.0001)
    FROM test600 b
    WHERE st_intersects(a.geom, b.geom)
         and a.area > 0
         and b.area > 0
         AND a.area > b.area
         and wof_id is not null
         )
      ), 
      a.geom) 
    FROM test600 a 
    where a.area > 0 
    ;