xref: /netbsd-src/external/bsd/ntp/dist/util/kern.c (revision cdfa2a7ef92791ba9db70a584a1d904730e6fb46)
1 /*	$NetBSD: kern.c,v 1.5 2020/05/25 20:47:37 christos Exp $	*/
2 
3 /*
4  * This program simulates a first-order, type-II phase-lock loop using
5  * actual code segments from modified kernel distributions for SunOS,
6  * Ultrix and OSF/1 kernels. These segments do not use any licensed code.
7  */
8 
9 #ifdef HAVE_CONFIG_H
10 # include <config.h>
11 #endif
12 
13 #include <stdio.h>
14 #include <ctype.h>
15 #include <math.h>
16 #include <sys/time.h>
17 
18 #ifdef HAVE_TIMEX_H
19 # include "timex.h"
20 #endif
21 
22 /*
23  * Phase-lock loop definitions
24  */
25 #define HZ 100			/* timer interrupt frequency (Hz) */
26 #define MAXPHASE 512000		/* max phase error (us) */
27 #define MAXFREQ 200		/* max frequency error (ppm) */
28 #define TAU 2			/* time constant (shift 0 - 6) */
29 #define POLL 16			/* interval between updates (s) */
30 #define MAXSEC 1200		/* max interval between updates (s) */
31 
32 /*
33  * Function declarations
34  */
35 void hardupdate();
36 void hardclock();
37 void second_overflow();
38 
39 /*
40  * Kernel variables
41  */
42 int tick;			/* timer interrupt period (us) */
43 int fixtick;			/* amortization constant (ppm) */
44 struct timeval timex;		/* ripoff of kernel time variable */
45 
46 /*
47  * Phase-lock loop variables
48  */
49 int time_status = TIME_BAD;	/* clock synchronization status */
50 long time_offset = 0;		/* time adjustment (us) */
51 long time_constant = 0;		/* pll time constant */
52 long time_tolerance = MAXFREQ;	/* frequency tolerance (ppm) */
53 long time_precision = 1000000 / HZ; /* clock precision (us) */
54 long time_maxerror = MAXPHASE;	/* maximum error (us) */
55 long time_esterror = MAXPHASE;	/* estimated error (us) */
56 long time_phase = 0;		/* phase offset (scaled us) */
57 long time_freq = 0;		/* frequency offset (scaled ppm) */
58 long time_adj = 0;		/* tick adjust (scaled 1 / HZ) */
59 long time_reftime = 0;		/* time at last adjustment (s) */
60 
61 /*
62  * Simulation variables
63  */
64 double timey = 0;		/* simulation time (us) */
65 long timez = 0;			/* current error (us) */
66 long poll_interval = 0;		/* poll counter */
67 
68 /*
69  * Simulation test program
70  */
71 int
main(int argc,char * argv[])72 main(
73 	int argc,
74 	char *argv[]
75 	)
76 {
77 	tick = 1000000 / HZ;
78 	fixtick = 1000000 % HZ;
79 	timex.tv_sec = 0;
80 	timex.tv_usec = MAXPHASE;
81 	time_freq = 0;
82 	time_constant = TAU;
83 	printf("tick %d us, fixtick %d us\n", tick, fixtick);
84 	printf("      time    offset      freq   _offset     _freq      _adj\n");
85 
86 	/*
87 	 * Grind the loop until ^C
88 	 */
89 	while (1) {
90 		timey += (double)(1000000) / HZ;
91 		if (timey >= 1000000)
92 		    timey -= 1000000;
93 		hardclock();
94 		if (timex.tv_usec >= 1000000) {
95 			timex.tv_usec -= 1000000;
96 			timex.tv_sec++;
97 			second_overflow();
98 			poll_interval++;
99 			if (!(poll_interval % POLL)) {
100 				timez = (long)timey - timex.tv_usec;
101 				if (timez > 500000)
102 				    timez -= 1000000;
103 				if (timez < -500000)
104 				    timez += 1000000;
105 				hardupdate(timez);
106 				printf("%10li%10li%10.2f  %08lx  %08lx  %08lx\n",
107 				       timex.tv_sec, timez,
108 				       (double)time_freq / (1 << SHIFT_KF),
109 				       time_offset, time_freq, time_adj);
110 			}
111 		}
112 	}
113 }
114 
115 /*
116  * This routine simulates the ntp_adjtime() call
117  *
118  * For default SHIFT_UPDATE = 12, offset is limited to +-512 ms, the
119  * maximum interval between updates is 4096 s and the maximum frequency
120  * offset is +-31.25 ms/s.
121  */
122 void
hardupdate(long offset)123 hardupdate(
124 	long offset
125 	)
126 {
127 	long ltemp, mtemp;
128 
129 	time_offset = offset << SHIFT_UPDATE;
130 	mtemp = timex.tv_sec - time_reftime;
131 	time_reftime = timex.tv_sec;
132 	if (mtemp > MAXSEC)
133 	    mtemp = 0;
134 
135 	/* ugly multiply should be replaced */
136 	if (offset < 0)
137 	    time_freq -= (-offset * mtemp) >>
138 		    (time_constant + time_constant);
139 	else
140 	    time_freq += (offset * mtemp) >>
141 		    (time_constant + time_constant);
142 	ltemp = time_tolerance << SHIFT_KF;
143 	if (time_freq > ltemp)
144 	    time_freq = ltemp;
145 	else if (time_freq < -ltemp)
146 	    time_freq = -ltemp;
147 	if (time_status == TIME_BAD)
148 	    time_status = TIME_OK;
149 }
150 
151 /*
152  * This routine simulates the timer interrupt
153  */
154 void
hardclock(void)155 hardclock(void)
156 {
157 	int ltemp, time_update;
158 
159 	time_update = tick;	/* computed by adjtime() */
160 	time_phase += time_adj;
161 	if (time_phase < -FINEUSEC) {
162 		ltemp = -time_phase >> SHIFT_SCALE;
163 		time_phase += ltemp << SHIFT_SCALE;
164 		time_update -= ltemp;
165 	}
166 	else if (time_phase > FINEUSEC) {
167 		ltemp = time_phase >> SHIFT_SCALE;
168 		time_phase -= ltemp << SHIFT_SCALE;
169 		time_update += ltemp;
170 	}
171 	timex.tv_usec += time_update;
172 }
173 
174 /*
175  * This routine simulates the overflow of the microsecond field
176  *
177  * With SHIFT_SCALE = 23, the maximum frequency adjustment is +-256 us
178  * per tick, or 25.6 ms/s at a clock frequency of 100 Hz. The time
179  * contribution is shifted right a minimum of two bits, while the frequency
180  * contribution is a right shift. Thus, overflow is prevented if the
181  * frequency contribution is limited to half the maximum or 15.625 ms/s.
182  */
183 void
second_overflow(void)184 second_overflow(void)
185 {
186 	int ltemp;
187 
188 	time_maxerror += time_tolerance;
189 	if (time_offset < 0) {
190 		ltemp = -time_offset >>
191 			(SHIFT_KG + time_constant);
192 		time_offset += ltemp;
193 		time_adj = -(ltemp <<
194 			     (SHIFT_SCALE - SHIFT_HZ - SHIFT_UPDATE));
195 	} else {
196 		ltemp = time_offset >>
197 			(SHIFT_KG + time_constant);
198 		time_offset -= ltemp;
199 		time_adj = ltemp <<
200 			(SHIFT_SCALE - SHIFT_HZ - SHIFT_UPDATE);
201 	}
202 	if (time_freq < 0)
203 	    time_adj -= -time_freq >> (SHIFT_KF + SHIFT_HZ - SHIFT_SCALE);
204 	else
205 	    time_adj += time_freq >> (SHIFT_KF + SHIFT_HZ - SHIFT_SCALE);
206 	time_adj += fixtick << (SHIFT_SCALE - SHIFT_HZ);
207 
208 	/* ugly divide should be replaced */
209 	if (timex.tv_sec % 86400 == 0) {
210 		switch (time_status) {
211 
212 		    case TIME_INS:
213 			timex.tv_sec--; /* !! */
214 			time_status = TIME_OOP;
215 			break;
216 
217 		    case TIME_DEL:
218 			timex.tv_sec++;
219 			time_status = TIME_OK;
220 			break;
221 
222 		    case TIME_OOP:
223 			time_status = TIME_OK;
224 			break;
225 		}
226 	}
227 }
228