Skip site navigation (1)Skip section navigation (2)
Date:      Sun, 12 Aug 2012 23:10:32 -0000
From:      Bruce Evans <brde@optusnet.com.au>
To:        Stephen Montgomery-Smith <stephen@missouri.edu>
Cc:        Diane Bruce <db@db.net>, Bruce Evans <brde@optusnet.com.au>, John Baldwin <jhb@freebsd.org>, David Chisnall <theraven@freebsd.org>, Bruce Evans <bde@freebsd.org>, Steve Kargl <sgk@troutmask.apl.washington.edu>, David Schultz <das@freebsd.org>, Peter Jeremy <peter@rulingia.com>, Warner Losh <imp@bsdimp.com>
Subject:   Re: Use of C99 extra long double math functions after r236148
Message-ID:  <20120723044308.X6145@besplex.bde.org>
Resent-Message-ID: <20120812231025.GX20453@server.rulingia.com>
In-Reply-To: <50095CDE.4050507@missouri.edu>
References:  <20120714120432.GA70706@server.rulingia.com> <20120717084457.U3890@besplex.bde.org> <5004A5C7.1040405@missouri.edu> <5004DEA9.1050001@missouri.edu> <20120717200931.U6624@besplex.bde.org> <5006D13D.2080702@missouri.edu> <20120718205625.GA409@troutmask.apl.washington.edu> <500725F2.7060603@missouri.edu> <20120719025345.GA1376@troutmask.apl.washington.edu> <50077987.1080307@missouri.edu> <20120719032706.GA1558@troutmask.apl.washington.edu> <5007826D.7060806@missouri.edu> <5007AD41.9070000@missouri.edu> <20120719205347.T2601@besplex.bde.org> <50084322.7020401@missouri.edu> <20120720035001.W4053@besplex.bde.org> <50085441.4090305@missouri.edu> <20120720162953.N2162@besplex.bde.org> <20120720184114.B2790@besplex.bde.org> <50095CDE.4050507@missouri.edu>

next in thread | previous in thread | raw e-mail | index | archive | help
  This message is in MIME format.  The first part should be readable text,
  while the remaining parts are likely unreadable without MIME-aware tools.

--0-1467877092-1342985376=:6145
Content-Type: TEXT/PLAIN; charset=US-ASCII; format=flowed

Replying again to this...

On Fri, 20 Jul 2012, Stephen Montgomery-Smith wrote:

> I worked quite a bit on my clog function last night.  Could you look at this 
> one?

I ended up deleting most of your changes in this one.

> The trick for when hypot(x,y) is close to 1 can only work so well, and you 
> are testing it out of its range of applicability.  But for the special case x 
> is close to 1, and y is close to zero, I have a much better formula.  I can 
> produce three other formula so that I can handle |x| close to 1, y small, and 
> |y| close to 1 and x small.
>
> After fixing it up, could you send it back as an attachment?  That will make 
> it easier for me to put it back into my system, and work more on it.

It will be painful for everyone to understand and merge.  This time as
an attachment as well as (partly) inline with commentary.

% #include <complex.h>
% #include <float.h>
% #include <math.h>
% 
% #include "math_private.h"
% 
% static inline void
% xonorm(double *a, double *b)
% {
% 	double w;
% 
% 	if (fabs(*a) < fabs(*b)) {
% 		w = *a;
% 		*a = *b;
% 		*b = w;
% 	}
% 	STRICT_ASSIGN(double, w, *a + *b);
% 	*b = (*a - w) + *b;
% 	*a = w;
% }
% 
% #define	xnorm(a, b)	xonorm(&(a), &(b))
% 
% #define	xspadd(a, b, c) do {	\
% 	double __tmp;		\
% 				\
% 	__tmp = (c);		\
% 	xonorm(&__tmp, &(a));	\
% 	(b) += (a);		\
% 	(a) = __tmp;		\
% } while (0)
% 
% static inline void
% xonormf(float *a, float *b)
% {
% 	float w;
% 
% 	if (fabsf(*a) < fabsf(*b)) {
% 		w = *a;
% 		*a = *b;
% 		*b = w;
% 	}
% 	STRICT_ASSIGN(float, w, *a + *b);
% 	*b = (*a - w) + *b;
% 	*a = w;
% }
% 
% #define	xnormf(a, b)	xonormf(&(a), &(b))
% 
% #define	xspaddf(a, b, c) do {	\
% 	float __tmp;		\
% 				\
% 	__tmp = (c);		\
% 	xonormf(&__tmp, &(a));	\
% 	(b) += (a);		\
% 	(a) = __tmp;		\
% } while (0)

Above are my standard extra-precision macros from my math_private.h,
cut down for use here and named with an x.  Then expanded and pessimized
to swap the args.  Optimal callers ensure that they don't need swapping.
See s_fma.c for a fuller algorithm that doesn't need swapping but does
more operations.  I started to copy that, but s_fma.c doesn't seem to
have anything as convenient as xspaddf().  Further expanded and
pessimized() to do STRICT_ASSIGN().  Optimal callers use float_t and
double_t so that STRICT_ASSIGN() is unnecessary.  Compiler bugs break
algorithms like the above on i386 unless float_t and double_t, or
STRICT_ASSIGN() are used.  Fixing the bugs would give the same slowness
as STRICT_ASSIGN(), but globally by doing it for all assignments, so
even I now consider these bugs to be features and C standards to be
broken for not allowing them.  I first discussed fixing them with gcc
maintainers over 20 years ago.

% 
% double complex
% clog(double complex z)
% {
% 	double x, y, h, t1, t2, t3;
% 	double ax, ay, x0, y0, x1, y1;
% 
% 	x = creal(z);
% 	y = cimag(z);
% 
% 	/* Handle NaNs using the general formula to mix them right. */
% 	if (x != x || y != y)
% 		return (cpack(log(hypot(x, y)), atan2(y, x)));

I replaced all my messy ifdefs for this by the function call.  Also
changes isnan() to a not-so-magic test.  Though isnan() is about the
only FP classification macro that I trust the compiler for.

% 
% 	ax = fabs(x);
% 	ay = fabs(y);
% 	if (ax < ay) {
% 		t1 = ax;
% 		ax = ay;
% 		ay = t1;
% 	}

I got tired of repeating fabs()'s, and need to know which arg is larger.

% 
% 	/*
% 	 * To avoid unnecessary overflow, if x or y are very large, divide x
% 	 * and y by M_E, and then add 1 to the logarithm.  This depends on
% 	 * M_E being larger than sqrt(2).
% 	 *
% 	 * XXX bugs to fix:
% 	 * - underflow if one of x or y is tiny.  e_hypot.c avoids this
% 	 *   problem, and optimizes for the case that the ratio of the
% 	 *   args is very large, by returning the absolute value of
% 	 *   the largest arg in this case.
% 	 * - not very accurate.  Could divide by 2 and add log(2) in extra
% 	 *   precision.  A general scaling step that divides by 2**k and
% 	 *   adds k*log(2) in extra precision might be good for reducing
% 	 *   the range so that we don't have to worry about overflow or
% 	 *   underflow in the general steps.  This needs the previous step
% 	 *   of eliminating large ratios of args so that the args can be
% 	 *   scaled on the same scale.
% 	 * - s/are/is/ in comment.
% 	 */
% 	if (ax > 1e308)
% 		return (cpack(log(hypot(x / M_E, y / M_E)) + 1, atan2(y, x)));

No need to avoid infinities here.  No need to test y now that we know
ax is largest.

% 
% 	if (ax == 1) {
% 		if (ay < 1e-150)
% 			return (cpack((ay * 0.5) * ay, atan2(y, x)));
% 		return (cpack(log1p(ay * ay) * 0.5, atan2(y, x)));
% 	}

Special case mainly for when (ay * ay) is rounded down to the smallest
denormal.

% 
% 	/*
% 	 * Because atan2 and hypot conform to C99, this also covers all the
% 	 * edge cases when x or y are 0 or infinite.
% 	 */
% 	if (ax < 1e-50 || ay < 1e-50 || ax > 1e50 || ay > 1e50)
% 		return (cpack(log(hypot(x, y)), atan2(y, x)));

Not quite right.  It is the ratio that matters more than these magic
magnitudes.

% 
% 	/* We don't need to worry about overflow in x*x+y*y. */
% 
% 	/*
% 	 * Take extra care so that ULP of real part is small if h is
% 	 * moderately close to 1.  If one only cares about the relative error
% 	 * of the whole result (real and imaginary part taken together), this
% 	 * algorithm is overkill.
% 	 *
% 	 * This algorithm does a rather good job if |h-1| >= 1e-5.  The only
% 	 * algorithm that I can think of that would work for any h close to
% 	 * one would require hypot(x,y) being computed using double double
% 	 * precision precision (i.e. double as many bits in the mantissa as
% 	 * double precision).
% 	 *
% 	 * x0 and y0 are good approximations to x and y, but have their bits
% 	 * trimmed so that double precision floating point is capable of
% 	 * calculating x0*x0 + y0*y0 - 1 exactly.
% 	 */

Comments not all updated.  This one especially out of date.

% 	x0 = ax;
% 	SET_LOW_WORD(x0, 0);
% 	x1 = ax - x0;
% 	y0 = ay;
% 	SET_LOW_WORD(y0, 0);
% 	y1 = ay - y0;

Sloppy decomposition with only 21 bits in hi part.  Since we are short
of bits, we shouldn't burn 5 like this for efficency.  In float precision,
all the multiplications are exact since 24 splits exactly.

% 	/* Notice that mathematically, h = t1*(1+t3). */
% #if 0

Old version.  Still drops bits unnecessary, although I added several
full hi/lo decomposition steps to it.

% 	t1 = x0 * x0;		/* Exact. */
% 	t2 = y0 * y0;		/* Exact. */

Comments not quite right.  All of the muliplications are as exact as
possible.  They would all be exact if we could split 53 in half, and
did so.

% 	STRICT_ASSIGN(double, t3, t1 + t2);
% 	t2 = (t1 - t3) + t2;
% 	t1 = t3;		/* Now t1+t2 is hi+lo for x0*x0+y0*y0.*/
% 	t2 += 2 * x0 * x1;
% 	STRICT_ASSIGN(double, t3, t1 + t2);
% 	t2 = (t1 - t3) + t2;
% 	t1 = t3;
% 	t2 += 2 * y0 * y1;
% 	STRICT_ASSIGN(double, t3, t1 + t2);
% 	t2 = (t1 - t3) + t2;
% 	t1 = t3;
% 	t2 += x1 * x1 + y1 * y1;
% 	STRICT_ASSIGN(double, t3, t1 + t2);
% 	t2 = (t1 - t3) + t2;
% 	t1 = t3;		/* Now t1+t2 is hi+lo for x*x+y*y.*/
% #else
% 	t1 = x1 * x1;
% 	t2 = y1 * y1;
% 	xnorm(t1, t2);
% 	xspadd(t1, t2, 2 * y0 * y1);
% 	xspadd(t1, t2, 2 * x0 * x1);
% 	xspadd(t1, t2, y0 * y0);
% 	xspadd(t1, t2, x0 * x0);
% 	xnorm(t1, t2);

It was too hard to turn the above into this without using the macros.
Now all the multiplications are as exact as possible, and extra precision
is used for all the additions (this mattered even for the first 2 terms).
Terms should be added from the smallest to the highest.  This happens in
most cases and some bits are lost when it isn't.

% #endif
% 	t3 = t2 / t1;
% 	/*
% 	 * |t3| ~< 2**-22 since we work with 24 extra bits of precision, so
% 	 * log1p(t3) can be evaluated with about 13 extra bits of precision
% 	 * using 2 terms of its power series.  But there are complexities
% 	 * to avoid underflow.
% 	 */

Complexities to avoid underflow incomplete and not here yet.

Comment otherwise anachronistic/anaspac(sp?)istic.  22 and 13 are for the
float version.  The final xnorm() step (to maximize accuracy) ensures that
|t3| < 2**-24 precisely for floats (half an ulp).  2**-53 for doubles.

% 	return (cpack((t3 - t3*0.5*t3 + log(t1)) * 0.5, atan2(y, x)));

The second term for log1p(t3) is probably nonsense .  We lose by
inaccuracies in t3 itself.

% }
% 
% float complex
% clogf(float complex z)
% {

This is a routine translation.  It duplicates too many comments.  Only
this has been tested much with the latest accuracty fixes.

% 	float x, y, h, t1, t2, t3;
% 	float ax, ay, x0, y0, x1, y1;
% 	uint32_t hx, hy;
% 
% 	x = crealf(z);
% 	y = cimagf(z);
% 
% 	/* Handle NaNs using the general formula to mix them right. */
% 	if (x != x || y != y)
% 		return (cpack(log(hypot(x, y)), atan2(y, x)));

Oops, copied too much -- forgot to add f's.

Bruce
--0-1467877092-1342985376=:6145
Content-Type: TEXT/PLAIN; charset=US-ASCII; name="cplex.c"
Content-Transfer-Encoding: BASE64
Content-ID: <20120723052936.D6145@besplex.bde.org>
Content-Description: 
Content-Disposition: attachment; filename="cplex.c"

I2luY2x1ZGUgPGNvbXBsZXguaD4NCiNpbmNsdWRlIDxmbG9hdC5oPg0KI2lu
Y2x1ZGUgPG1hdGguaD4NCg0KI2luY2x1ZGUgIm1hdGhfcHJpdmF0ZS5oIg0K
DQpzdGF0aWMgaW5saW5lIHZvaWQNCnhvbm9ybShkb3VibGUgKmEsIGRvdWJs
ZSAqYikNCnsNCglkb3VibGUgdzsNCg0KCWlmIChmYWJzKCphKSA8IGZhYnMo
KmIpKSB7DQoJCXcgPSAqYTsNCgkJKmEgPSAqYjsNCgkJKmIgPSB3Ow0KCX0N
CglTVFJJQ1RfQVNTSUdOKGRvdWJsZSwgdywgKmEgKyAqYik7DQoJKmIgPSAo
KmEgLSB3KSArICpiOw0KCSphID0gdzsNCn0NCg0KI2RlZmluZQl4bm9ybShh
LCBiKQl4b25vcm0oJihhKSwgJihiKSkNCg0KI2RlZmluZQl4c3BhZGQoYSwg
YiwgYykgZG8gewlcDQoJZG91YmxlIF9fdG1wOwkJXA0KCQkJCVwNCglfX3Rt
cCA9IChjKTsJCVwNCgl4b25vcm0oJl9fdG1wLCAmKGEpKTsJXA0KCShiKSAr
PSAoYSk7CQlcDQoJKGEpID0gX190bXA7CQlcDQp9IHdoaWxlICgwKQ0KDQpz
dGF0aWMgaW5saW5lIHZvaWQNCnhvbm9ybWYoZmxvYXQgKmEsIGZsb2F0ICpi
KQ0Kew0KCWZsb2F0IHc7DQoNCglpZiAoZmFic2YoKmEpIDwgZmFic2YoKmIp
KSB7DQoJCXcgPSAqYTsNCgkJKmEgPSAqYjsNCgkJKmIgPSB3Ow0KCX0NCglT
VFJJQ1RfQVNTSUdOKGZsb2F0LCB3LCAqYSArICpiKTsNCgkqYiA9ICgqYSAt
IHcpICsgKmI7DQoJKmEgPSB3Ow0KfQ0KDQojZGVmaW5lCXhub3JtZihhLCBi
KQl4b25vcm1mKCYoYSksICYoYikpDQoNCiNkZWZpbmUJeHNwYWRkZihhLCBi
LCBjKSBkbyB7CVwNCglmbG9hdCBfX3RtcDsJCVwNCgkJCQlcDQoJX190bXAg
PSAoYyk7CQlcDQoJeG9ub3JtZigmX190bXAsICYoYSkpOwlcDQoJKGIpICs9
IChhKTsJCVwNCgkoYSkgPSBfX3RtcDsJCVwNCn0gd2hpbGUgKDApDQoNCmRv
dWJsZSBjb21wbGV4DQpjbG9nKGRvdWJsZSBjb21wbGV4IHopDQp7DQoJZG91
YmxlIHgsIHksIGgsIHQxLCB0MiwgdDM7DQoJZG91YmxlIGF4LCBheSwgeDAs
IHkwLCB4MSwgeTE7DQoNCgl4ID0gY3JlYWwoeik7DQoJeSA9IGNpbWFnKHop
Ow0KDQoJLyogSGFuZGxlIE5hTnMgdXNpbmcgdGhlIGdlbmVyYWwgZm9ybXVs
YSB0byBtaXggdGhlbSByaWdodC4gKi8NCglpZiAoeCAhPSB4IHx8IHkgIT0g
eSkNCgkJcmV0dXJuIChjcGFjayhsb2coaHlwb3QoeCwgeSkpLCBhdGFuMih5
LCB4KSkpOw0KDQoJYXggPSBmYWJzKHgpOw0KCWF5ID0gZmFicyh5KTsNCglp
ZiAoYXggPCBheSkgew0KCQl0MSA9IGF4Ow0KCQlheCA9IGF5Ow0KCQlheSA9
IHQxOw0KCX0NCg0KCS8qDQoJICogVG8gYXZvaWQgdW5uZWNlc3Nhcnkgb3Zl
cmZsb3csIGlmIHggb3IgeSBhcmUgdmVyeSBsYXJnZSwgZGl2aWRlIHgNCgkg
KiBhbmQgeSBieSBNX0UsIGFuZCB0aGVuIGFkZCAxIHRvIHRoZSBsb2dhcml0
aG0uICBUaGlzIGRlcGVuZHMgb24NCgkgKiBNX0UgYmVpbmcgbGFyZ2VyIHRo
YW4gc3FydCgyKS4NCgkgKg0KCSAqIFhYWCBidWdzIHRvIGZpeDoNCgkgKiAt
IHVuZGVyZmxvdyBpZiBvbmUgb2YgeCBvciB5IGlzIHRpbnkuICBlX2h5cG90
LmMgYXZvaWRzIHRoaXMNCgkgKiAgIHByb2JsZW0sIGFuZCBvcHRpbWl6ZXMg
Zm9yIHRoZSBjYXNlIHRoYXQgdGhlIHJhdGlvIG9mIHRoZQ0KCSAqICAgYXJn
cyBpcyB2ZXJ5IGxhcmdlLCBieSByZXR1cm5pbmcgdGhlIGFic29sdXRlIHZh
bHVlIG9mDQoJICogICB0aGUgbGFyZ2VzdCBhcmcgaW4gdGhpcyBjYXNlLg0K
CSAqIC0gbm90IHZlcnkgYWNjdXJhdGUuICBDb3VsZCBkaXZpZGUgYnkgMiBh
bmQgYWRkIGxvZygyKSBpbiBleHRyYQ0KCSAqICAgcHJlY2lzaW9uLiAgQSBn
ZW5lcmFsIHNjYWxpbmcgc3RlcCB0aGF0IGRpdmlkZXMgYnkgMioqayBhbmQN
CgkgKiAgIGFkZHMgaypsb2coMikgaW4gZXh0cmEgcHJlY2lzaW9uIG1pZ2h0
IGJlIGdvb2QgZm9yIHJlZHVjaW5nDQoJICogICB0aGUgcmFuZ2Ugc28gdGhh
dCB3ZSBkb24ndCBoYXZlIHRvIHdvcnJ5IGFib3V0IG92ZXJmbG93IG9yDQoJ
ICogICB1bmRlcmZsb3cgaW4gdGhlIGdlbmVyYWwgc3RlcHMuICBUaGlzIG5l
ZWRzIHRoZSBwcmV2aW91cyBzdGVwDQoJICogICBvZiBlbGltaW5hdGluZyBs
YXJnZSByYXRpb3Mgb2YgYXJncyBzbyB0aGF0IHRoZSBhcmdzIGNhbiBiZQ0K
CSAqICAgc2NhbGVkIG9uIHRoZSBzYW1lIHNjYWxlLg0KCSAqIC0gcy9hcmUv
aXMvIGluIGNvbW1lbnQuDQoJICovDQoJaWYgKGF4ID4gMWUzMDgpDQoJCXJl
dHVybiAoY3BhY2sobG9nKGh5cG90KHggLyBNX0UsIHkgLyBNX0UpKSArIDEs
IGF0YW4yKHksIHgpKSk7DQoNCglpZiAoYXggPT0gMSkgew0KCQlpZiAoYXkg
PCAxZS0xNTApDQoJCQlyZXR1cm4gKGNwYWNrKChheSAqIDAuNSkgKiBheSwg
YXRhbjIoeSwgeCkpKTsNCgkJcmV0dXJuIChjcGFjayhsb2cxcChheSAqIGF5
KSAqIDAuNSwgYXRhbjIoeSwgeCkpKTsNCgl9DQoNCgkvKg0KCSAqIEJlY2F1
c2UgYXRhbjIgYW5kIGh5cG90IGNvbmZvcm0gdG8gQzk5LCB0aGlzIGFsc28g
Y292ZXJzIGFsbCB0aGUNCgkgKiBlZGdlIGNhc2VzIHdoZW4geCBvciB5IGFy
ZSAwIG9yIGluZmluaXRlLg0KCSAqLw0KCWlmIChheCA8IDFlLTUwIHx8IGF5
IDwgMWUtNTAgfHwgYXggPiAxZTUwIHx8IGF5ID4gMWU1MCkNCgkJcmV0dXJu
IChjcGFjayhsb2coaHlwb3QoeCwgeSkpLCBhdGFuMih5LCB4KSkpOw0KDQoJ
LyogV2UgZG9uJ3QgbmVlZCB0byB3b3JyeSBhYm91dCBvdmVyZmxvdyBpbiB4
KngreSp5LiAqLw0KDQoJLyoNCgkgKiBUYWtlIGV4dHJhIGNhcmUgc28gdGhh
dCBVTFAgb2YgcmVhbCBwYXJ0IGlzIHNtYWxsIGlmIGggaXMNCgkgKiBtb2Rl
cmF0ZWx5IGNsb3NlIHRvIDEuICBJZiBvbmUgb25seSBjYXJlcyBhYm91dCB0
aGUgcmVsYXRpdmUgZXJyb3INCgkgKiBvZiB0aGUgd2hvbGUgcmVzdWx0IChy
ZWFsIGFuZCBpbWFnaW5hcnkgcGFydCB0YWtlbiB0b2dldGhlciksIHRoaXMN
CgkgKiBhbGdvcml0aG0gaXMgb3ZlcmtpbGwuDQoJICoNCgkgKiBUaGlzIGFs
Z29yaXRobSBkb2VzIGEgcmF0aGVyIGdvb2Qgam9iIGlmIHxoLTF8ID49IDFl
LTUuICBUaGUgb25seQ0KCSAqIGFsZ29yaXRobSB0aGF0IEkgY2FuIHRoaW5r
IG9mIHRoYXQgd291bGQgd29yayBmb3IgYW55IGggY2xvc2UgdG8NCgkgKiBv
bmUgd291bGQgcmVxdWlyZSBoeXBvdCh4LHkpIGJlaW5nIGNvbXB1dGVkIHVz
aW5nIGRvdWJsZSBkb3VibGUNCgkgKiBwcmVjaXNpb24gcHJlY2lzaW9uIChp
LmUuIGRvdWJsZSBhcyBtYW55IGJpdHMgaW4gdGhlIG1hbnRpc3NhIGFzDQoJ
ICogZG91YmxlIHByZWNpc2lvbikuDQoJICoNCgkgKiB4MCBhbmQgeTAgYXJl
IGdvb2QgYXBwcm94aW1hdGlvbnMgdG8geCBhbmQgeSwgYnV0IGhhdmUgdGhl
aXIgYml0cw0KCSAqIHRyaW1tZWQgc28gdGhhdCBkb3VibGUgcHJlY2lzaW9u
IGZsb2F0aW5nIHBvaW50IGlzIGNhcGFibGUgb2YNCgkgKiBjYWxjdWxhdGlu
ZyB4MCp4MCArIHkwKnkwIC0gMSBleGFjdGx5Lg0KCSAqLw0KCXgwID0gYXg7
DQoJU0VUX0xPV19XT1JEKHgwLCAwKTsNCgl4MSA9IGF4IC0geDA7DQoJeTAg
PSBheTsNCglTRVRfTE9XX1dPUkQoeTAsIDApOw0KCXkxID0gYXkgLSB5MDsN
CgkvKiBOb3RpY2UgdGhhdCBtYXRoZW1hdGljYWxseSwgaCA9IHQxKigxK3Qz
KS4gKi8NCiNpZiAwDQoJdDEgPSB4MCAqIHgwOwkJLyogRXhhY3QuICovDQoJ
dDIgPSB5MCAqIHkwOwkJLyogRXhhY3QuICovDQoJU1RSSUNUX0FTU0lHTihk
b3VibGUsIHQzLCB0MSArIHQyKTsNCgl0MiA9ICh0MSAtIHQzKSArIHQyOw0K
CXQxID0gdDM7CQkvKiBOb3cgdDErdDIgaXMgaGkrbG8gZm9yIHgwKngwK3kw
KnkwLiovDQoJdDIgKz0gMiAqIHgwICogeDE7DQoJU1RSSUNUX0FTU0lHTihk
b3VibGUsIHQzLCB0MSArIHQyKTsNCgl0MiA9ICh0MSAtIHQzKSArIHQyOw0K
CXQxID0gdDM7DQoJdDIgKz0gMiAqIHkwICogeTE7DQoJU1RSSUNUX0FTU0lH
Tihkb3VibGUsIHQzLCB0MSArIHQyKTsNCgl0MiA9ICh0MSAtIHQzKSArIHQy
Ow0KCXQxID0gdDM7DQoJdDIgKz0geDEgKiB4MSArIHkxICogeTE7DQoJU1RS
SUNUX0FTU0lHTihkb3VibGUsIHQzLCB0MSArIHQyKTsNCgl0MiA9ICh0MSAt
IHQzKSArIHQyOw0KCXQxID0gdDM7CQkvKiBOb3cgdDErdDIgaXMgaGkrbG8g
Zm9yIHgqeCt5KnkuKi8NCiNlbHNlDQoJdDEgPSB4MSAqIHgxOw0KCXQyID0g
eTEgKiB5MTsNCgl4bm9ybSh0MSwgdDIpOw0KCXhzcGFkZCh0MSwgdDIsIDIg
KiB5MCAqIHkxKTsNCgl4c3BhZGQodDEsIHQyLCAyICogeDAgKiB4MSk7DQoJ
eHNwYWRkKHQxLCB0MiwgeTAgKiB5MCk7DQoJeHNwYWRkKHQxLCB0MiwgeDAg
KiB4MCk7DQoJeG5vcm0odDEsIHQyKTsNCiNlbmRpZg0KCXQzID0gdDIgLyB0
MTsNCgkvKg0KCSAqIHx0M3wgfjwgMioqLTIyIHNpbmNlIHdlIHdvcmsgd2l0
aCAyNCBleHRyYSBiaXRzIG9mIHByZWNpc2lvbiwgc28NCgkgKiBsb2cxcCh0
MykgY2FuIGJlIGV2YWx1YXRlZCB3aXRoIGFib3V0IDEzIGV4dHJhIGJpdHMg
b2YgcHJlY2lzaW9uDQoJICogdXNpbmcgMiB0ZXJtcyBvZiBpdHMgcG93ZXIg
c2VyaWVzLiAgQnV0IHRoZXJlIGFyZSBjb21wbGV4aXRpZXMNCgkgKiB0byBh
dm9pZCB1bmRlcmZsb3cuDQoJICovDQoJcmV0dXJuIChjcGFjaygodDMgLSB0
MyowLjUqdDMgKyBsb2codDEpKSAqIDAuNSwgYXRhbjIoeSwgeCkpKTsNCn0N
Cg0KZmxvYXQgY29tcGxleA0KY2xvZ2YoZmxvYXQgY29tcGxleCB6KQ0Kew0K
CWZsb2F0IHgsIHksIGgsIHQxLCB0MiwgdDM7DQoJZmxvYXQgYXgsIGF5LCB4
MCwgeTAsIHgxLCB5MTsNCgl1aW50MzJfdCBoeCwgaHk7DQoNCgl4ID0gY3Jl
YWxmKHopOw0KCXkgPSBjaW1hZ2Yoeik7DQoNCgkvKiBIYW5kbGUgTmFOcyB1
c2luZyB0aGUgZ2VuZXJhbCBmb3JtdWxhIHRvIG1peCB0aGVtIHJpZ2h0LiAq
Lw0KCWlmICh4ICE9IHggfHwgeSAhPSB5KQ0KCQlyZXR1cm4gKGNwYWNrKGxv
ZyhoeXBvdCh4LCB5KSksIGF0YW4yKHksIHgpKSk7DQoNCglheCA9IGZhYnNm
KHgpOw0KCWF5ID0gZmFic2YoeSk7DQoJaWYgKGF4IDwgYXkpIHsNCgkJdDEg
PSBheDsNCgkJYXggPSBheTsNCgkJYXkgPSB0MTsNCgl9DQoNCgkvKg0KCSAq
IFRvIGF2b2lkIHVubmVjZXNzYXJ5IG92ZXJmbG93LCBpZiB4IG9yIHkgYXJl
IHZlcnkgbGFyZ2UsIGRpdmlkZSB4DQoJICogYW5kIHkgYnkgTV9FLCBhbmQg
dGhlbiBhZGQgMSB0byB0aGUgbG9nYXJpdGhtLiAgVGhpcyBkZXBlbmRzIG9u
DQoJICogTV9FIGJlaW5nIGxhcmdlciB0aGFuIHNxcnQoMikuDQoJICoNCgkg
KiBYWFggYnVncyB0byBmaXg6DQoJICogLSB1bmRlcmZsb3cgaWYgb25lIG9m
IHggb3IgeSBpcyB0aW55LiAgZV9oeXBvdC5jIGF2b2lkcyB0aGlzDQoJICog
ICBwcm9ibGVtLCBhbmQgb3B0aW1pemVzIGZvciB0aGUgY2FzZSB0aGF0IHRo
ZSByYXRpbyBvZiB0aGUNCgkgKiAgIGFyZ3MgaXMgdmVyeSBsYXJnZSwgYnkg
cmV0dXJuaW5nIHRoZSBhYnNvbHV0ZSB2YWx1ZSBvZg0KCSAqICAgdGhlIGxh
cmdlc3QgYXJnIGluIHRoaXMgY2FzZS4NCgkgKiAtIG5vdCB2ZXJ5IGFjY3Vy
YXRlLiAgQ291bGQgZGl2aWRlIGJ5IDIgYW5kIGFkZCBsb2coMikgaW4gZXh0
cmENCgkgKiAgIHByZWNpc2lvbi4gIEEgZ2VuZXJhbCBzY2FsaW5nIHN0ZXAg
dGhhdCBkaXZpZGVzIGJ5IDIqKmsgYW5kDQoJICogICBhZGRzIGsqbG9nKDIp
IGluIGV4dHJhIHByZWNpc2lvbiBtaWdodCBiZSBnb29kIGZvciByZWR1Y2lu
Zw0KCSAqICAgdGhlIHJhbmdlIHNvIHRoYXQgd2UgZG9uJ3QgaGF2ZSB0byB3
b3JyeSBhYm91dCBvdmVyZmxvdyBvcg0KCSAqICAgdW5kZXJmbG93IGluIHRo
ZSBnZW5lcmFsIHN0ZXBzLiAgVGhpcyBuZWVkcyB0aGUgcHJldmlvdXMgc3Rl
cA0KCSAqICAgb2YgZWxpbWluYXRpbmcgbGFyZ2UgcmF0aW9zIG9mIGFyZ3Mg
c28gdGhhdCB0aGUgYXJncyBjYW4gYmUNCgkgKiAgIHNjYWxlZCBvbiB0aGUg
c2FtZSBzY2FsZS4NCgkgKiAtIHMvYXJlL2lzLyBpbiBjb21tZW50Lg0KCSAq
Lw0KCWlmIChheCA+IDFlMzhGKQ0KCQlyZXR1cm4gKGNwYWNrKGxvZ2YoaHlw
b3RmKHggLyAoZmxvYXQpTV9FLCB5IC8gKGZsb2F0KU1fRSkpICsgMSwNCgkJ
ICAgIGF0YW4yZih5LCB4KSkpOw0KDQoJaWYgKGF4ID09IDEpIHsNCgkJaWYg
KGF5IDwgMWUtMThGKQ0KCQkJcmV0dXJuIChjcGFja2YoKGF5ICogMC41Rikg
KiBheSwgYXRhbjJmKHksIHgpKSk7DQoJCXJldHVybiAoY3BhY2tmKGxvZzFw
ZihheSAqIGF5KSAqIDAuNUYsIGF0YW4yZih5LCB4KSkpOw0KCX0NCg0KCS8q
DQoJICogQmVjYXVzZSBhdGFuMiBhbmQgaHlwb3QgY29uZm9ybSB0byBDOTks
IHRoaXMgYWxzbyBjb3ZlcnMgYWxsIHRoZQ0KCSAqIGVkZ2UgY2FzZXMgd2hl
biB4IG9yIHkgYXJlIDAgb3IgaW5maW5pdGUuDQoJICovDQoJaWYgKGF4IDwg
MWUtMTBGIHx8IGF5IDwgMWUtMTBGIHx8IGF4ID4gMWUxMEYgfHwgYXkgPiAx
ZTEwRikNCgkJcmV0dXJuIChjcGFja2YobG9nZihoeXBvdGYoeCwgeSkpLCBh
dGFuMmYoeSwgeCkpKTsNCg0KCS8qIFdlIGRvbid0IG5lZWQgdG8gd29ycnkg
YWJvdXQgb3ZlcmZsb3cgaW4geCp4K3kqeS4gKi8NCg0KCS8qDQoJICogVGFr
ZSBleHRyYSBjYXJlIHNvIHRoYXQgVUxQIG9mIHJlYWwgcGFydCBpcyBzbWFs
bCBpZiBoIGlzDQoJICogbW9kZXJhdGVseSBjbG9zZSB0byAxLiAgSWYgb25l
IG9ubHkgY2FyZXMgYWJvdXQgdGhlIHJlbGF0aXZlIGVycm9yDQoJICogb2Yg
dGhlIHdob2xlIHJlc3VsdCAocmVhbCBhbmQgaW1hZ2luYXJ5IHBhcnQgdGFr
ZW4gdG9nZXRoZXIpLCB0aGlzDQoJICogYWxnb3JpdGhtIGlzIG92ZXJraWxs
Lg0KCSAqDQoJICogVGhpcyBhbGdvcml0aG0gZG9lcyBhIHJhdGhlciBnb29k
IGpvYiBpZiB8aC0xfCA+PSAxZS01LiAgVGhlIG9ubHkNCgkgKiBhbGdvcml0
aG0gdGhhdCBJIGNhbiB0aGluayBvZiB0aGF0IHdvdWxkIHdvcmsgZm9yIGFu
eSBoIGNsb3NlIHRvDQoJICogb25lIHdvdWxkIHJlcXVpcmUgaHlwb3QoeCx5
KSBiZWluZyBjb21wdXRlZCB1c2luZyBkb3VibGUgZG91YmxlDQoJICogcHJl
Y2lzaW9uIHByZWNpc2lvbiAoaS5lLiBkb3VibGUgYXMgbWFueSBiaXRzIGlu
IHRoZSBtYW50aXNzYSBhcw0KCSAqIGRvdWJsZSBwcmVjaXNpb24pLg0KCSAq
DQoJICogeDAgYW5kIHkwIGFyZSBnb29kIGFwcHJveGltYXRpb25zIHRvIHgg
YW5kIHksIGJ1dCBoYXZlIHRoZWlyIGJpdHMNCgkgKiB0cmltbWVkIHNvIHRo
YXQgZG91YmxlIHByZWNpc2lvbiBmbG9hdGluZyBwb2ludCBpcyBjYXBhYmxl
IG9mDQoJICogY2FsY3VsYXRpbmcgeDAqeDAgKyB5MCp5MCAtIDEgZXhhY3Rs
eS4NCgkgKi8NCg0KCUdFVF9GTE9BVF9XT1JEKGh4LCB4KTsNCglTRVRfRkxP
QVRfV09SRCh4MCwgaHggJiAweGZmZmZmMDAwKTsNCgl4MSA9IHggLSB4MDsN
CglHRVRfRkxPQVRfV09SRChoeSwgeSk7DQoJU0VUX0ZMT0FUX1dPUkQoeTAs
IGh5ICYgMHhmZmZmZjAwMCk7DQoJeTEgPSB5IC0geTA7DQoJLyogTm90aWNl
IHRoYXQgbWF0aGVtYXRpY2FsbHksIGggPSB0MSooMSt0MykuICovDQojaWYg
MA0KCXQxID0geDAgKiB4MDsJCS8qIEV4YWN0LiAqLw0KCXQyID0geTAgKiB5
MDsJCS8qIEV4YWN0LiAqLw0KCVNUUklDVF9BU1NJR04oZmxvYXQsIHQzLCB0
MSArIHQyKTsNCgl0MiA9ICh0MSAtIHQzKSArIHQyOw0KCXQxID0gdDM7CQkv
KiBOb3cgdDErdDIgaXMgaGkrbG8gZm9yIHgwKngwK3kwKnkwLiovDQoJdDIg
Kz0gMiAqIHgwICogeDE7DQoJU1RSSUNUX0FTU0lHTihmbG9hdCwgdDMsIHQx
ICsgdDIpOw0KCXQyID0gKHQxIC0gdDMpICsgdDI7DQoJdDEgPSB0MzsNCgl0
MiArPSAyICogeTAgKiB5MTsNCglTVFJJQ1RfQVNTSUdOKGZsb2F0LCB0Mywg
dDEgKyB0Mik7DQoJdDIgPSAodDEgLSB0MykgKyB0MjsNCgl0MSA9IHQzOw0K
CXQyICs9IHgxICogeDEgKyB5MSAqIHkxOw0KCVNUUklDVF9BU1NJR04oZmxv
YXQsIHQzLCB0MSArIHQyKTsNCgl0MiA9ICh0MSAtIHQzKSArIHQyOw0KCXQx
ID0gdDM7CQkvKiBOb3cgdDErdDIgaXMgaGkrbG8gZm9yIHgqeCt5KnkuKi8N
CiNlbHNlDQoJdDEgPSB4MSAqIHgxOw0KCXQyID0geTEgKiB5MTsNCgl4bm9y
bWYodDEsIHQyKTsNCgl4c3BhZGRmKHQxLCB0MiwgMiAqIHkwICogeTEpOw0K
CXhzcGFkZGYodDEsIHQyLCAyICogeDAgKiB4MSk7DQoJeHNwYWRkZih0MSwg
dDIsIHkwICogeTApOw0KCXhzcGFkZGYodDEsIHQyLCB4MCAqIHgwKTsNCgl4
bm9ybWYodDEsIHQyKTsNCiNlbmRpZg0KCXQzID0gdDIgLyB0MTsNCgkvKg0K
CSAqIHx0M3wgfjwgMioqLTEwIHNpbmNlIHdlIHdvcmsgd2l0aCAxMiBleHRy
YSBiaXRzIG9mIHByZWNpc2lvbiwgc28NCgkgKiBsb2cxcCh0MykgY2FuIGJl
IGV2YWx1YXRlZCB3aXRoIGFib3V0IDYgZXh0cmEgYml0cyBvZiBwcmVjaXNp
b24NCgkgKiB1c2luZyAyIHRlcm1zIG9mIGl0cyBwb3dlciBzZXJpZXMuICBC
dXQgdGhlcmUgYXJlIGNvbXBsZXhpdGllcw0KCSAqIHRvIGF2b2lkIHVuZGVy
Zmxvdy4NCgkgKi8NCglyZXR1cm4gKGNwYWNrZigodDMgLSB0MyowLjVGKnQz
ICsgbG9nZih0MSkpICogMC41RiwgYXRhbjJmKHksIHgpKSk7DQp9DQo=

--0-1467877092-1342985376=:6145--



Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20120723044308.X6145>