Well, I found the bug. There was an initialization error which made the
program think there was a shorter path than there was in a number of cases.
Alas, as always, the bug did not reveal itself in the test runs I did, it
was only in the final totals that it was apparent.

The bad thing is that removing the bug increases the compute time considerably,
because you have in essence to calculate an explicit path for every
configuration. I estimate the increase is by a factor of about 3, based on
some experiments.

The good thing is that it got me thinking about an old idea I had.
To calculate path lengths to the group generated by [U,D,L2,R2,F2,B2]
it is irrelevant whether you rotate the complete cube around the U-D
axis. Also half turn rotations around the R-L axis are irrelevant.
And finally, mirroring the cube along the F-R-B-L plane is irrelevant!
But in most cases this changes twist/flip and position values. So if
we look at the twist/flip/position space of 2187*2048*495 cosets we can
reduce the calculations along one dimension as long as we remember the
number of representations. I did some calculations and found that
2187 can be reduced to 168 or 2048 can be reduced to 219 or 495 can be
reduced to 45. Of these obviously the first one is the most fruitful.
Of course I have to check myself (and anybody who is willing, go ahead).
Reduction of 2187 to 168 would reduce the total space from 2,217,093,120
to 170,311,680 calculations. But then again, dividing the latter figure
by 4, that could be done in core on a machine with about 50 MByte of memory
(to calculate path lengths in core only 2 bits per configuration are needed).

I will try around, dik