# [metapost] Turning Number and QuickSort

luigi scarso luigi.scarso at gmail.com
Wed Jan 21 09:14:27 CET 2015

```It seems that ￼the University College Cork has some problem to deliver
emails to the metapost list,
so I forward some notes that
Marc van Dongen has sent to me last days:

"""
I have noticed lots of problems with standard metapost, including
lots of problems related to rounding and an error in turningnumber.

I re-implemented turningnumber (in mpost) and it seems to work for
my purpose. The implementation is quite simple and it avoids rounding
errors. (It's still possible to reduce more rounding but I don't
think it's needed.) Please let me know if you want the implementation.
(I still have to do some proper testing on it.)
"""
The files are in attachment.

--
luigi
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://tug.org/pipermail/metapost/attachments/20150121/562008e0/attachment.html>
-------------- next part --------------
% This is a simple implementation to determine whether a  path
% is drawn anti clockwise. To do this, it computes the signed
% area of the path consisting of the support points of the path.
% ASSUMPTION: _pat_ does have a non-zero area.
vardef IsAntiClockwise( expr pat ) =
((cycle pat) and (Area( pat ) > AREA_ACCURACY))
%(turningnumber( pat ) > 0)
enddef;

vardef Area( expr pat ) =
save a, len, sum; numeric a[], len[], sum;
len := length pat;
% First we compute (twice) the signed areas of the triangles
% _(0,0) -- (point i of pat) -- (point (i + 1) of pat) -- cycle_
for i = len - 1 downto 0:
a[ i ] := Det( point i of pat, point (i + 1) of pat );
endfor
% Since adding the sizes may result in large rounding errors,
% we sort them in descending order
SortAbsGeq( a )( len );
% The remaining statements compute the sum of the signed areas
% in a ``robust'' fashion, trying to minimise rounding. A little
% bit more robust would be to first compute the frequencies of
% the entries in _a_, then multiplying each unique number in
% _a_ by its frequency, then sorting the products, and finally
% e.g. [1,2,1,1];
% sort it -> [2,1,1,1];
% compute frequencies -> [2:1,1:3];
% compute products -> [2,3];
% sort -> [2,3];
sum := 0;
for i = len - 1 downto 0:
sum := sum + a[ i ];
endfor
% We could have just returned _sum_ but
% this way we get the actual signed area.
0.5sum
enddef;

endinput
-------------- next part --------------
vardef QuickSort( suffix a, cmp )( expr size ) =
QuickSort_( a, cmp )( 0, size - 1 );
enddef;

vardef QuickSort_( suffix a, cmp )( expr lo, hi ) =
if lo < hi:
save sep; numeric sep;
%Swop( a )( lo, floor( (lo + hi) / 2 ) );
Swop( a )( lo, lo + floor( uniformdeviate( hi - lo - 1 ) ) );
sep := Partition( a, cmp, lo, hi );
Swop( a )( lo, sep );
QuickSort_( a, cmp )( lo, sep - 1 );
QuickSort_( a, cmp )( sep + 1, hi );
fi
enddef;

vardef Partition( suffix a, cmp )( expr lo, hi ) =
save dest, pivot; numeric dest, pivot;

dest := lo;
pivot := a[ lo ];
for i = lo + 1 upto hi:
if cmp( a[ i ], pivot ):
dest := dest + 1;
Swop( a )( dest, i );
fi
endfor
dest
enddef;

vardef Swop( suffix a )( expr fst, snd ) =
save tmp; numeric tmp;
tmp := a[ fst ];
a[ fst ] := a[ snd ];
a[ snd ] := tmp;
enddef;

vardef SortAbsGeq( suffix a )( expr size ) =
save CMP__;
vardef CMP__( expr fst, snd ) =
(abs( fst ) >= abs( snd ))
enddef;
QuickSort( a, CMP__ )( size );
enddef;

endinput
```