Home | History | Annotate | Line # | Download | only in timed
      1 /*	$NetBSD: networkdelta.c,v 1.13 2021/10/18 14:16:49 christos Exp $	*/
      2 
      3 /*-
      4  * Copyright (c) 1985, 1993 The Regents of the University of California.
      5  * All rights reserved.
      6  *
      7  * Redistribution and use in source and binary forms, with or without
      8  * modification, are permitted provided that the following conditions
      9  * are met:
     10  * 1. Redistributions of source code must retain the above copyright
     11  *    notice, this list of conditions and the following disclaimer.
     12  * 2. Redistributions in binary form must reproduce the above copyright
     13  *    notice, this list of conditions and the following disclaimer in the
     14  *    documentation and/or other materials provided with the distribution.
     15  * 3. Neither the name of the University nor the names of its contributors
     16  *    may be used to endorse or promote products derived from this software
     17  *    without specific prior written permission.
     18  *
     19  * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
     20  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
     21  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
     22  * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
     23  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
     24  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
     25  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
     26  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
     27  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
     28  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
     29  * SUCH DAMAGE.
     30  */
     31 
     32 #include <sys/cdefs.h>
     33 #ifndef lint
     34 #if 0
     35 static char sccsid[] = "@(#)networkdelta.c	8.3 (Berkeley) 4/27/95";
     36 #else
     37 __RCSID("$NetBSD: networkdelta.c,v 1.13 2021/10/18 14:16:49 christos Exp $");
     38 #endif
     39 #endif /* not lint */
     40 
     41 #include "globals.h"
     42 
     43 static long median(float, float*, long*, long*, unsigned int);
     44 
     45 /*
     46  * Compute a corrected date.
     47  *	Compute the median of the reasonable differences.  First compute
     48  *	the median of all authorized differences, and then compute the
     49  *	median of all differences that are reasonably close to the first
     50  *	median.
     51  *
     52  * This differs from the original BSD implementation, which looked for
     53  *	the largest group of machines with essentially the same date.
     54  *	That assumed that machines with bad clocks would be uniformly
     55  *	distributed.  Unfortunately, in real life networks, the distribution
     56  *	of machines is not uniform among models of machines, and the
     57  *	distribution of errors in clocks tends to be quite consistent
     58  *	for a given model.  In other words, all model VI Supre Servres
     59  *	from GoFast Inc. tend to have about the same error.
     60  *	The original BSD implementation would chose the clock of the
     61  *	most common model, and discard all others.
     62  *
     63  *	Therefore, get best we can do is to try to average over all
     64  *	of the machines in the network, while discarding "obviously"
     65  *	bad values.
     66  */
     67 long
     68 networkdelta(void)
     69 {
     70 	struct hosttbl *htp;
     71 	long med;
     72 	long lodelta, hidelta;
     73 	long logood, higood;
     74 	long x[NHOSTS];
     75 	long *xp;
     76 	int numdelta;
     77 	float eps;
     78 
     79 	/*
     80 	 * compute the median of the good values
     81 	 */
     82 	med = 0;
     83 	numdelta = 1;
     84 	xp = &x[0];
     85 	*xp = 0;			/* account for ourself */
     86 	for (htp = self.l_fwd; htp != &self; htp = htp->l_fwd) {
     87 		if (htp->good
     88 		    && htp->noanswer == 0
     89 		    && htp->delta != HOSTDOWN) {
     90 			med += htp->delta;
     91 			numdelta++;
     92 			*++xp = htp->delta;
     93 		}
     94 	}
     95 
     96 	/*
     97 	 * If we are the only trusted time keeper, then do not change our
     98 	 * clock.  There may be another time keeping service active.
     99 	 */
    100 	if (numdelta == 1)
    101 		return 0;
    102 
    103 	med /= numdelta;
    104 	eps = med - x[0];
    105 	if (trace)
    106 		fprintf(fd, "median of %d values starting at %ld is about ",
    107 			numdelta, med);
    108 	med = median(med, &eps, &x[0], xp+1, VALID_RANGE);
    109 
    110 	/*
    111 	 * compute the median of all values near the good median
    112 	 */
    113 	hidelta = med + GOOD_RANGE;
    114 	lodelta = med - GOOD_RANGE;
    115 	higood = med + VGOOD_RANGE;
    116 	logood = med - VGOOD_RANGE;
    117 	xp = &x[0];
    118 	htp = &self;
    119 	do {
    120 		if (htp->noanswer == 0
    121 		    && htp->delta >= lodelta
    122 		    && htp->delta <= hidelta
    123 		    && (htp->good
    124 			|| (htp->delta >= logood
    125 			    && htp->delta <= higood))) {
    126 			*xp++ = htp->delta;
    127 		}
    128 	} while (&self != (htp = htp->l_fwd));
    129 
    130 	if (xp == &x[0]) {
    131 		if (trace)
    132 			fprintf(fd, "nothing close to median %ld\n", med);
    133 		return med;
    134 	}
    135 
    136 	if (xp == &x[1]) {
    137 		if (trace)
    138 			fprintf(fd, "only value near median is %ld\n", x[0]);
    139 		return x[0];
    140 	}
    141 
    142 	if (trace)
    143 		fprintf(fd, "median of %ld values starting at %ld is ",
    144 		        (long)(xp - &x[0]), med);
    145 	return median(med, &eps, &x[0], xp, 1);
    146 }
    147 
    148 
    149 /*
    150  * compute the median of an array of signed integers, using the idea
    151  *	in <<Numerical Recipes>>.
    152  */
    153 static long
    154 median(float a,				/* initial guess for the median */
    155        float *eps_ptr,			/* spacing near the median */
    156        long *x, long *xlim,		/* the data */
    157        unsigned int gnuf)		/* good enough estimate */
    158 {
    159 	long *xptr;
    160 	float ap = (float)LONG_MAX;	/* bounds on the median */
    161 	float am = (float)-LONG_MAX;
    162 	float aa;
    163 	int npts;			/* # of points above & below guess */
    164 	float xp;			/* closet point above the guess */
    165 	float xm;			/* closet point below the guess */
    166 	float eps;
    167 	float dum, sum, sumx;
    168 	int pass;
    169 #define AMP	1.5			/* smoothing constants */
    170 #define AFAC	1.5
    171 
    172 	eps = *eps_ptr;
    173 	if (eps < 1.0) {
    174 		eps = -eps;
    175 		if (eps < 1.0)
    176 			eps = 1.0;
    177 	}
    178 
    179 	for (pass = 1; ; pass++) {	/* loop over the data */
    180 		sum = 0.0;
    181 		sumx = 0.0;
    182 		npts = 0;
    183 		xp = (float)LONG_MAX;
    184 		xm = (float)-LONG_MAX;
    185 
    186 		for (xptr = x; xptr != xlim; xptr++) {
    187 			float xx = *xptr;
    188 
    189 			dum = xx - a;
    190 			if (dum != 0.0) {	/* avoid dividing by 0 */
    191 				if (dum > 0.0) {
    192 					npts++;
    193 					if (xx < xp)
    194 						xp = xx;
    195 				} else {
    196 					npts--;
    197 					if (xx > xm)
    198 						xm = xx;
    199 					dum = -dum;
    200 				}
    201 				dum = 1.0/(eps + dum);
    202 				sum += dum;
    203 				sumx += xx * dum;
    204 			}
    205 		}
    206 
    207 		if (ap-am < gnuf || sum == 0) {
    208 			if (trace)
    209 				fprintf(fd,
    210 			           "%ld in %d passes; early out balance=%d\n",
    211 				        (long)a, pass, npts);
    212 			return a;	/* guess was good enough */
    213 		}
    214 
    215 		aa = (sumx/sum-a)*AMP;
    216 		if (npts >= 2) {	/* guess was too low */
    217 			am = a;
    218 			aa = xp + max(0.0, aa);
    219 			if (aa > ap)
    220 				aa = (a + ap)/2;
    221 
    222 		} else if (npts <= -2) {  /* guess was two high */
    223 			ap = a;
    224 			aa = xm + min(0.0, aa);
    225 			if (aa < am)
    226 				aa = (a + am)/2;
    227 
    228 		} else {
    229 			break;		/* got it */
    230 		}
    231 
    232 		if (a == aa) {
    233 			if (trace)
    234 				fprintf(fd,
    235 				  "%ld in %d passes; force out balance=%d\n",
    236 				        (long)a, pass, npts);
    237 			return a;
    238 		}
    239 		eps = AFAC*abs(aa - a);
    240 		*eps_ptr = eps;
    241 		a = aa;
    242 	}
    243 
    244 	if (((x - xlim) % 2) != 0) {    /* even number of points? */
    245 		if (npts == 0)		/* yes, return an average */
    246 			a = (xp+xm)/2;
    247 		else if (npts > 0)
    248 			a =  (a+xp)/2;
    249 		else
    250 			a = (xm+a)/2;
    251 
    252 	} else 	if (npts != 0) {	/* odd number of points */
    253 		if (npts > 0)
    254 			a = xp;
    255 		else
    256 			a = xm;
    257 	}
    258 
    259 	if (trace)
    260 		fprintf(fd, "%ld in %d passes\n", (long)a, pass);
    261 	return a;
    262 }
    263