/******************************************************************************/
/* */
/* COMPUTE THE MINIMUM ELEMENT IN A LOOP */
/* 04/15/10 (dkc) */
/* */
/* This C program finds the minimum element in a loop. c is set to 1 or -1. */
/* Potential loops are in limbs in S of a least-residue tree. n is the */
/* number of odd elements in the limb. The periods for divisors (d) of n */
/* are; */
/* */
/* d=5, p=13 (or 14) */
/* d=7, p=10 (or 11) */
/* d=12, p=51 (or 52) */
/* d=17, p=18 (or 17) (large period, 17=2*5+7) */
/* d=19, p=8 (or 9) */
/* d=22, p=8 (or 7) */
/* d=26, p=5 (or 4) */
/* d=27, p=5 (or 4) */
/* d=29, p=28 (or 27) (large period, 29=3*5+2*7) */
/* d=31, p=8 (or 7) */
/* d=32, p=3 (or 4) */
/* d=33, p=3 (or 4) */
/* d=37, p=3 (or 2) */
/* d=39, p=5 (or 6) */
/* d=41, p=60 (or 61) (large period, 41=4*5+3*7) */
/* d=43, p=6 (or 7) */
/* d=46, p=10 (or 11) */
/* d=47, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=53, p=331 (or ?) (large period, 53=5*5+4*7) */
/* d=59, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=61, p=3 (or 4) */
/* d=67, p=5 (or 6) */
/* d=69, p=3 (or 2) */
/* d=71, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=73, p=3 (or 4) */
/* d=79, p=5 (or 4) */
/* d=81, p=3 (or 2) */
/* d=83, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=88, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=89, p=16 (or 17) (large period, 89=8*5+7*7) */
/* */
/* Periods for divisors that are multiples of the above divisors are; */
/* */
/* d=10, p=6 (or 7) */
/* d=14, p=5 (or 6) */
/* d=15, p=5 (or 4) */
/* d=20, p=3 (or 4) */
/* d=21, p=3 (or 4) */
/* d=24, p=26 (or 25) */
/* d=25, p=3 (or 2) */
/* d=28, p=3 (or 2) */
/* d=30, p=3 (or 2) */
/* d=34, p=8 (or 9) */
/* d=35, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=36, p=18 (or 17) */
/* d=38, p=4 (or 5) */
/* d=40, p=3 (or 2) */
/* d=42, p=3 (or 2) */
/* d=44, p=4 (or 3) */
/* d=45, p=3 (or 4) */
/* d=48, p=13 (or 12) */
/* d=49, p=3 (or 2) */
/* d=50, p=4 (or 5) */
/* d=51, p=5 (or 6) */
/* d=52, p=3 (or 2) */
/* d=54, p=3 (or 2) */
/* d=55, p=5 (or 6) */
/* d=56, p=4 (or 5) */
/* d=57, p=3 (or 2) */
/* d=58, p=13 (or 14) */
/* d=60, p=11 (or 10) */
/* d=62, p=4 (or 3) */
/* d=63, p=6 (or 7) */
/* d=64, p=3 (or 2) */
/* d=65, p=44 (or 45) (large period, 65=6*5+5*7) */
/* d=66, p=3 (or 2) */
/* d=68, p=4 (or 5) */
/* d=70, p=19 (or ?) */
/* d=72, p=9 (or 8) */
/* d=74, p=3 (or 4) */
/* d=75, p=7 (or 8) */
/* d=76, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=77, p=23 (or 24) (large period, 77=7*5+6*7) */
/* d=78, p=3 (or 2) */
/* d=80, p=4 (or 5) */
/* d=81, p=3 (or 2) */
/* d=82, p=30 (or 31) */
/* d=84, p=7 (or 8) */
/* d=85, p=3 (or 4) */
/* d=86, p=3 (or 4) */
/* d=87, p=9 (or 10) */
/* d=88, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=90, p=3 (or 2) */
/* */
/* Periods where d is not of the form 5x+7y are; */
/* */
/* d=2, p=6 (or 5) */
/* d=3, p=4 (or 5) */
/* d=4, p=3 (or 2) */
/* d=6, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=8, p=3 (or 4) */
/* d=9, p=4 (or 3) */
/* d=11, p=3 (or 2) */
/* d=13, p=3 (or 2) */
/* d=16, p=3 (or 2) */
/* d=18, p=2 (or 3) (sometimes oscillates above or below zero) */
/* d=23, p=3 (or 2) */
/* */
/* When the predominant period is greater than 2, the x values for a given */
/* d value have the general shape of a saw-tooth curve with spikes from */
/* either negative values to positive values or positive values to negative */
/* values. The longest periods (and hence the biggest spikes) occur when d */
/* is of the form 5*k+7*(k-1), k=1, 2, 3, .... */
/* */
/* The x values for prime pairs are of interest. Prime pairs (other than */
/* (3, 5)) can be classified into two types; in one type, the first prime */
/* of the pair is of the form 12*k+5 and in the other type, the first prime */
/* of the pair is of the form 12*k-1. The truncated x values for n=17, 29, */
/* 29, 41, 101, 137, 149, 197, 269, 281, 461, 521, 569, 617, 641, 809, 821, */
/* 857, and 881 are 108, 281, 867, -433, -340, -325, -288, -263, -261, 481, */
/* 770, 1263, 2740, 5699, -1196, -1121, -950, and -868 respectively. The */
/* truncated x values for n=19, 31, 43, 103, 139, 151, 199, 271, 283, 463, */
/* 523, 571, 619, 643, 811, 823, 859, and 883 are -58, -81, -98, -142, -156, */
/* -159, -169, -179, -180, 989, 2809, -14586, -2340, -1705, -695, -673, -617,*/
/* and -587 respectively. The truncated x values for n=11, 59, 71, 107, 179,*/
/* 191, 227, 239, 311, 347, 419, 431, 599, 659, and 827 are -9, -40, -46, */
/* 90, 213, 243, 369, 427, 1408, 6689, -1465, -1259, -529, -466, and 1215 */
/* respectively. The truncated x values for n=13, 61, 73, 109, 181, 193, */
/* 229, 241, 313, 349, 421, 433, 601, 661, and 829 are 11, 66, 85, 158, 516, */
/* 655, 1825, 3482, -1165, -797, -543, -520, -371, 674, and 4353 respect- */
/* ively. */
/* */
/******************************************************************************/
#include <stdio.h>
#include <math.h>
void shiftn(unsigned int *a, unsigned int *b, unsigned int shift, unsigned int n);
void lshiftn(unsigned int *a, unsigned int *b, unsigned int shift, unsigned int n);
void div128_64(unsigned int a0, unsigned int a1, unsigned int a2,
unsigned int a3, unsigned int *quotient, unsigned int d2,
unsigned int d3);
void addn(unsigned int *a, unsigned int *b, unsigned int n);
void subn(unsigned int *c, unsigned int *d, unsigned int n);
void copyn(unsigned int *a, unsigned int *b, unsigned int n);
void setn(unsigned int *a, unsigned int b, unsigned int n);
void shiftn(unsigned int *a, unsigned int *b, unsigned int shift, unsigned int n);
unsigned int orn(unsigned int *a, unsigned int n);
int d=5; // divisor of n (should be of the form 5x+7y)
unsigned int ab[2582*2]={ // a and b values
0, 1, 0, 2, 1, 2, 1, 3, 2, 3,
2, 4, 3, 4, 4, 4, 4, 5, 5, 5,
5, 6, 6, 6, 7, 6, 7, 7, 8, 7,
8, 8, 9, 8, 9, 9, 10, 9, 11, 9,
11, 10, 12, 10, 12, 11, 13, 11, 14, 11,
14, 12, 15, 12, 15, 13, 16, 13, 16, 14,
17, 14, 18, 14, 18, 15, 19, 15, 19, 16,
20, 16, 21, 16, 21, 17, 22, 17, 22, 18,
23, 18, 23, 19, 24, 19, 25, 19, 25, 20,
26, 20, 26, 21, 27, 21, 28, 21, 28, 22,
29, 22, 29, 23, 30, 23, 31, 23, 31, 24,
32, 24, 32, 25, 33, 25, 33, 26, 34, 26,
35, 26, 35, 27, 36, 27, 36, 28, 37, 28,
38, 28, 38, 29, 39, 29, 39, 30, 40, 30,
40, 31, 41, 31, 42, 31, 42, 32, 43, 32,
43, 33, 44, 33, 45, 33, 45, 34, 46, 34,
46, 35, 47, 35, 47, 36, 48, 36, 49, 36,
49, 37, 50, 37, 50, 38, 51, 38, 52, 38,
52, 39, 53, 39, 53, 40, 54, 40, 54, 41,
55, 41, 56, 41, 56, 42, 57, 42, 57, 43,
58, 43, 59, 43, 59, 44, 60, 44, 60, 45,
61, 45, 62, 45, 62, 46, 63, 46, 63, 47,
64, 47, 64, 48, 65, 48, 66, 48, 66, 49,
67, 49, 67, 50, 68, 50, 69, 50, 69, 51,
70, 51, 70, 52, 71, 52, 71, 53, 72, 53,
73, 53, 73, 54, 74, 54, 74, 55, 75, 55,
76, 55, 76, 56, 77, 56, 77, 57, 78, 57,
78, 58, 79, 58, 80, 58, 80, 59, 81, 59,
81, 60, 82, 60, 83, 60, 83, 61, 84, 61,
84, 62, 85, 62, 85, 63, 86, 63, 87, 63,
87, 64, 88, 64, 88, 65, 89, 65, 90, 65,
90, 66, 91, 66, 91, 67, 92, 67, 93, 67,
93, 68, 94, 68, 94, 69, 95, 69, 95, 70,
96, 70, 97, 70, 97, 71, 98, 71, 98, 72,
99, 72, 100, 72, 100, 73, 101, 73, 101, 74,
102, 74, 102, 75, 103, 75, 104, 75, 104, 76,
105, 76, 105, 77, 106, 77, 107, 77, 107, 78,
108, 78, 108, 79, 109, 79, 109, 80, 110, 80,
111, 80, 111, 81, 112, 81, 112, 82, 113, 82,
114, 82, 114, 83, 115, 83, 115, 84, 116, 84,
116, 85, 117, 85, 118, 85, 118, 86, 119, 86,
119, 87, 120, 87, 121, 87, 121, 88, 122, 88,
122, 89, 123, 89, 124, 89, 124, 90, 125, 90,
125, 91, 126, 91, 126, 92, 127, 92, 128, 92,
128, 93, 129, 93, 129, 94, 130, 94, 131, 94,
131, 95, 132, 95, 132, 96, 133, 96, 133, 97,
134, 97, 135, 97, 135, 98, 136, 98, 136, 99,
137, 99, 138, 99, 138, 100, 139, 100, 139, 101,
140, 101, 140, 102, 141, 102, 142, 102, 142, 103,
143, 103, 143, 104, 144, 104, 145, 104, 145, 105,
146, 105, 146, 106, 147, 106, 147, 107, 148, 107,
149, 107, 149, 108, 150, 108, 150, 109, 151, 109,
152, 109, 152, 110, 153, 110, 153, 111, 154, 111,
155, 111, 155, 112, 156, 112, 156, 113, 157, 113,
157, 114, 158, 114, 159, 114, 159, 115, 160, 115,
160, 116, 161, 116, 162, 116, 162, 117, 163, 117,
163, 118, 164, 118, 164, 119, 165, 119, 166, 119,
166, 120, 167, 120, 167, 121, 168, 121, 169, 121,
169, 122, 170, 122, 170, 123, 171, 123, 171, 124,
172, 124, 173, 124, 173, 125, 174, 125, 174, 126,
175, 126, 176, 126, 176, 127, 177, 127, 177, 128,
178, 128, 178, 129, 179, 129, 180, 129, 180, 130,
181, 130, 181, 131, 182, 131, 183, 131, 183, 132,
184, 132, 184, 133, 185, 133, 186, 133, 186, 134,
187, 134, 187, 135, 188, 135, 188, 136, 189, 136,
190, 136, 190, 137, 191, 137, 191, 138, 192, 138,
193, 138, 193, 139, 194, 139, 194, 140, 195, 140,
195, 141, 196, 141, 197, 141, 197, 142, 198, 142,
198, 143, 199, 143, 200, 143, 200, 144, 201, 144,
201, 145, 202, 145, 202, 146, 203, 146, 204, 146,
204, 147, 205, 147, 205, 148, 206, 148, 207, 148,
207, 149, 208, 149, 208, 150, 209, 150, 210, 150,
210, 151, 211, 151, 211, 152, 212, 152, 212, 153,
213, 153, 214, 153, 214, 154, 215, 154, 215, 155,
216, 155, 217, 155, 217, 156, 218, 156, 218, 157,
219, 157, 219, 158, 220, 158, 221, 158, 221, 159,
222, 159, 222, 160, 223, 160, 224, 160, 224, 161,
225, 161, 225, 162, 226, 162, 226, 163, 227, 163,
228, 163, 228, 164, 229, 164, 229, 165, 230, 165,
231, 165, 231, 166, 232, 166, 232, 167, 233, 167,
233, 168, 234, 168, 235, 168, 235, 169, 236, 169,
236, 170, 237, 170, 238, 170, 238, 171, 239, 171,
239, 172, 240, 172, 241, 172, 241, 173, 242, 173,
242, 174, 243, 174, 243, 175, 244, 175, 245, 175,
245, 176, 246, 176, 246, 177, 247, 177, 248, 177,
248, 178, 249, 178, 249, 179, 250, 179, 250, 180,
251, 180, 252, 180, 252, 181, 253, 181, 253, 182,
254, 182, 255, 182, 255, 183, 256, 183, 256, 184,
257, 184, 257, 185, 258, 185, 259, 185, 259, 186,
260, 186, 260, 187, 261, 187, 262, 187, 262, 188,
263, 188, 263, 189, 264, 189, 264, 190, 265, 190,
266, 190, 266, 191, 267, 191, 267, 192, 268, 192,
269, 192, 269, 193, 270, 193, 270, 194, 271, 194,
272, 194, 272, 195, 273, 195, 273, 196, 274, 196,
274, 197, 275, 197, 276, 197, 276, 198, 277, 198,
277, 199, 278, 199, 279, 199, 279, 200, 280, 200,
280, 201, 281, 201, 281, 202, 282, 202, 283, 202,
283, 203, 284, 203, 284, 204, 285, 204, 286, 204,
286, 205, 287, 205, 287, 206, 288, 206, 288, 207,
289, 207, 290, 207, 290, 208, 291, 208, 291, 209,
292, 209, 293, 209, 293, 210, 294, 210, 294, 211,
295, 211, 295, 212, 296, 212, 297, 212, 297, 213,
298, 213, 298, 214, 299, 214, 300, 214, 300, 215,
301, 215, 301, 216, 302, 216, 303, 216, 303, 217,
304, 217, 304, 218, 305, 218, 305, 219, 306, 219,
307, 219, 307, 220, 308, 220, 308, 221, 309, 221,
310, 221, 310, 222, 311, 222, 311, 223, 312, 223,
312, 224, 313, 224, 314, 224, 314, 225, 315, 225,
315, 226, 316, 226, 317, 226, 317, 227, 318, 227,
318, 228, 319, 228, 319, 229, 320, 229, 321, 229,
321, 230, 322, 230, 322, 231, 323, 231, 324, 231,
324, 232, 325, 232, 325, 233, 326, 233, 326, 234,
327, 234, 328, 234, 328, 235, 329, 235, 329, 236,
330, 236, 331, 236, 331, 237, 332, 237, 332, 238,
333, 238, 334, 238, 334, 239, 335, 239, 335, 240,
336, 240, 336, 241, 337, 241, 338, 241, 338, 242,
339, 242, 339, 243, 340, 243, 341, 243, 341, 244,
342, 244, 342, 245, 343, 245, 343, 246, 344, 246,
345, 246, 345, 247, 346, 247, 346, 248, 347, 248,
348, 248, 348, 249, 349, 249, 349, 250, 350, 250,
350, 251, 351, 251, 352, 251, 352, 252, 353, 252,
353, 253, 354, 253, 355, 253, 355, 254, 356, 254,
356, 255, 357, 255, 357, 256, 358, 256, 359, 256,
359, 257, 360, 257, 360, 258, 361, 258, 362, 258,
362, 259, 363, 259, 363, 260, 364, 260, 365, 260,
365, 261, 366, 261, 366, 262, 367, 262, 367, 263,
368, 263, 369, 263, 369, 264, 370, 264, 370, 265,
371, 265, 372, 265, 372, 266, 373, 266, 373, 267,
374, 267, 374, 268, 375, 268, 376, 268, 376, 269,
377, 269, 377, 270, 378, 270, 379, 270, 379, 271,
380, 271, 380, 272, 381, 272, 381, 273, 382, 273,
383, 273, 383, 274, 384, 274, 384, 275, 385, 275,
386, 275, 386, 276, 387, 276, 387, 277, 388, 277,
389, 277, 389, 278, 390, 278, 390, 279, 391, 279,
391, 280, 392, 280, 393, 280, 393, 281, 394, 281,
394, 282, 395, 282, 396, 282, 396, 283, 397, 283,
397, 284, 398, 284, 398, 285, 399, 285, 400, 285,
400, 286, 401, 286, 401, 287, 402, 287, 403, 287,
403, 288, 404, 288, 404, 289, 405, 289, 405, 290,
406, 290, 407, 290, 407, 291, 408, 291, 408, 292,
409, 292, 410, 292, 410, 293, 411, 293, 411, 294,
412, 294, 412, 295, 413, 295, 414, 295, 414, 296,
415, 296, 415, 297, 416, 297, 417, 297, 417, 298,
418, 298, 418, 299, 419, 299, 420, 299, 420, 300,
421, 300, 421, 301, 422, 301, 422, 302, 423, 302,
424, 302, 424, 303, 425, 303, 425, 304, 426, 304,
427, 304, 427, 305, 428, 305, 428, 306, 429, 306,
429, 307, 430, 307, 431, 307, 431, 308, 432, 308,
432, 309, 433, 309, 434, 309, 434, 310, 435, 310,
435, 311, 436, 311, 436, 312, 437, 312, 438, 312,
438, 313, 439, 313, 439, 314, 440, 314, 441, 314,
441, 315, 442, 315, 442, 316, 443, 316, 443, 317,
444, 317, 445, 317, 445, 318, 446, 318, 446, 319,
447, 319, 448, 319, 448, 320, 449, 320, 449, 321,
450, 321, 451, 321, 451, 322, 452, 322, 452, 323,
453, 323, 453, 324, 454, 324, 455, 324, 455, 325,
456, 325, 456, 326, 457, 326, 458, 326, 458, 327,
459, 327, 459, 328, 460, 328, 460, 329, 461, 329,
462, 329, 462, 330, 463, 330, 463, 331, 464, 331,
465, 331, 465, 332, 466, 332, 466, 333, 467, 333,
467, 334, 468, 334, 469, 334, 469, 335, 470, 335,
470, 336, 471, 336, 472, 336, 472, 337, 473, 337,
473, 338, 474, 338, 474, 339, 475, 339, 476, 339,
476, 340, 477, 340, 477, 341, 478, 341, 479, 341,
479, 342, 480, 342, 480, 343, 481, 343, 482, 343,
482, 344, 483, 344, 483, 345, 484, 345, 484, 346,
485, 346, 486, 346, 486, 347, 487, 347, 487, 348,
488, 348, 489, 348, 489, 349, 490, 349, 490, 350,
491, 350, 491, 351, 492, 351, 493, 351, 493, 352,
494, 352, 494, 353, 495, 353, 496, 353, 496, 354,
497, 354, 497, 355, 498, 355, 498, 356, 499, 356,
500, 356, 500, 357, 501, 357, 501, 358, 502, 358,
503, 358, 503, 359, 504, 359, 504, 360, 505, 360,
505, 361, 506, 361, 507, 361, 507, 362, 508, 362,
508, 363, 509, 363, 510, 363, 510, 364, 511, 364,
511, 365, 512, 365, 513, 365, 513, 366, 514, 366,
514, 367, 515, 367, 515, 368, 516, 368, 517, 368,
517, 369, 518, 369, 518, 370, 519, 370, 520, 370,
520, 371, 521, 371, 521, 372, 522, 372, 522, 373,
523, 373, 524, 373, 524, 374, 525, 374, 525, 375,
526, 375, 527, 375, 527, 376, 528, 376, 528, 377,
529, 377, 529, 378, 530, 378, 531, 378, 531, 379,
532, 379, 532, 380, 533, 380, 534, 380, 534, 381,
535, 381, 535, 382, 536, 382, 536, 383, 537, 383,
538, 383, 538, 384, 539, 384, 539, 385, 540, 385,
541, 385, 541, 386, 542, 386, 542, 387, 543, 387,
544, 387, 544, 388, 545, 388, 545, 389, 546, 389,
546, 390, 547, 390, 548, 390, 548, 391, 549, 391,
549, 392, 550, 392, 551, 392, 551, 393, 552, 393,
552, 394, 553, 394, 553, 395, 554, 395, 555, 395,
555, 396, 556, 396, 556, 397, 557, 397, 558, 397,
558, 398, 559, 398, 559, 399, 560, 399, 560, 400,
561, 400, 562, 400, 562, 401, 563, 401, 563, 402,
564, 402, 565, 402, 565, 403, 566, 403, 566, 404,
567, 404, 567, 405, 568, 405, 569, 405, 569, 406,
570, 406, 570, 407, 571, 407, 572, 407, 572, 408,
573, 408, 573, 409, 574, 409, 575, 409, 575, 410,
576, 410, 576, 411, 577, 411, 577, 412, 578, 412,
579, 412, 579, 413, 580, 413, 580, 414, 581, 414,
582, 414, 582, 415, 583, 415, 583, 416, 584, 416,
584, 417, 585, 417, 586, 417, 586, 418, 587, 418,
587, 419, 588, 419, 589, 419, 589, 420, 590, 420,
590, 421, 591, 421, 591, 422, 592, 422, 593, 422,
593, 423, 594, 423, 594, 424, 595, 424, 596, 424,
596, 425, 597, 425, 597, 426, 598, 426, 599, 426,
599, 427, 600, 427, 600, 428, 601, 428, 601, 429,
602, 429, 603, 429, 603, 430, 604, 430, 604, 431,
605, 431, 606, 431, 606, 432, 607, 432, 607, 433,
608, 433, 608, 434, 609, 434, 610, 434, 610, 435,
611, 435, 611, 436, 612, 436, 613, 436, 613, 437,
614, 437, 614, 438, 615, 438, 615, 439, 616, 439,
617, 439, 617, 440, 618, 440, 618, 441, 619, 441,
620, 441, 620, 442, 621, 442, 621, 443, 622, 443,
622, 444, 623, 444, 624, 444, 624, 445, 625, 445,
625, 446, 626, 446, 627, 446, 627, 447, 628, 447,
628, 448, 629, 448, 630, 448, 630, 449, 631, 449,
631, 450, 632, 450, 632, 451, 633, 451, 634, 451,
634, 452, 635, 452, 635, 453, 636, 453, 637, 453,
637, 454, 638, 454, 638, 455, 639, 455, 639, 456,
640, 456, 641, 456, 641, 457, 642, 457, 642, 458,
643, 458, 644, 458, 644, 459, 645, 459, 645, 460,
646, 460, 646, 461, 647, 461, 648, 461, 648, 462,
649, 462, 649, 463, 650, 463, 651, 463, 651, 464,
652, 464, 652, 465, 653, 465, 653, 466, 654, 466,
655, 466, 655, 467, 656, 467, 656, 468, 657, 468,
658, 468, 658, 469, 659, 469, 659, 470, 660, 470,
661, 470, 661, 471, 662, 471, 662, 472, 663, 472,
663, 473, 664, 473, 665, 473, 665, 474, 666, 474,
666, 475, 667, 475, 668, 475, 668, 476, 669, 476,
669, 477, 670, 477, 670, 478, 671, 478, 672, 478,
672, 479, 673, 479, 673, 480, 674, 480, 675, 480,
675, 481, 676, 481, 676, 482, 677, 482, 677, 483,
678, 483, 679, 483, 679, 484, 680, 484, 680, 485,
681, 485, 682, 485, 682, 486, 683, 486, 683, 487,
684, 487, 684, 488, 685, 488, 686, 488, 686, 489,
687, 489, 687, 490, 688, 490, 689, 490, 689, 491,
690, 491, 690, 492, 691, 492, 692, 492, 692, 493,
693, 493, 693, 494, 694, 494, 694, 495, 695, 495,
696, 495, 696, 496, 697, 496, 697, 497, 698, 497,
699, 497, 699, 498, 700, 498, 700, 499, 701, 499,
701, 500, 702, 500, 703, 500, 703, 501, 704, 501,
704, 502, 705, 502, 706, 502, 706, 503, 707, 503,
707, 504, 708, 504, 708, 505, 709, 505, 710, 505,
710, 506, 711, 506, 711, 507, 712, 507, 713, 507,
713, 508, 714, 508, 714, 509, 715, 509, 715, 510,
716, 510, 717, 510, 717, 511, 718, 511, 718, 512,
719, 512, 720, 512, 720, 513, 721, 513, 721, 514,
722, 514, 723, 514, 723, 515, 724, 515, 724, 516,
725, 516, 725, 517, 726, 517, 727, 517, 727, 518,
728, 518, 728, 519, 729, 519, 730, 519, 730, 520,
731, 520, 731, 521, 732, 521, 732, 522, 733, 522,
734, 522, 734, 523, 735, 523, 735, 524, 736, 524,
737, 524, 737, 525, 738, 525, 738, 526, 739, 526,
739, 527, 740, 527, 741, 527, 741, 528, 742, 528,
742, 529, 743, 529, 744, 529, 744, 530, 745, 530,
745, 531, 746, 531, 746, 532, 747, 532, 748, 532,
748, 533, 749, 533, 749, 534, 750, 534, 751, 534,
751, 535, 752, 535, 752, 536, 753, 536, 754, 536,
754, 537, 755, 537, 755, 538, 756, 538, 756, 539,
757, 539, 758, 539, 758, 540, 759, 540, 759, 541,
760, 541, 761, 541, 761, 542, 762, 542, 762, 543,
763, 543, 763, 544, 764, 544, 765, 544, 765, 545,
766, 545, 766, 546, 767, 546, 768, 546, 768, 547,
769, 547, 769, 548, 770, 548, 770, 549, 771, 549,
772, 549, 772, 550, 773, 550, 773, 551, 774, 551,
775, 551, 775, 552, 776, 552, 776, 553, 777, 553,
778, 553, 778, 554, 779, 554, 779, 555, 780, 555,
780, 556, 781, 556, 782, 556, 782, 557, 783, 557,
783, 558, 784, 558, 785, 558, 785, 559, 786, 559,
786, 560, 787, 560, 787, 561, 788, 561, 789, 561,
789, 562, 790, 562, 790, 563, 791, 563, 792, 563,
792, 564, 793, 564, 793, 565, 794, 565, 794, 566,
795, 566, 796, 566, 796, 567, 797, 567, 797, 568,
798, 568, 799, 568, 799, 569, 800, 569, 800, 570,
801, 570, 801, 571, 802, 571, 803, 571, 803, 572,
804, 572, 804, 573, 805, 573, 806, 573, 806, 574,
807, 574, 807, 575, 808, 575, 809, 575, 809, 576,
810, 576, 810, 577, 811, 577, 811, 578, 812, 578,
813, 578, 813, 579, 814, 579, 814, 580, 815, 580,
816, 580, 816, 581, 817, 581, 817, 582, 818, 582,
818, 583, 819, 583, 820, 583, 820, 584, 821, 584,
821, 585, 822, 585, 823, 585, 823, 586, 824, 586,
824, 587, 825, 587, 825, 588, 826, 588, 827, 588,
827, 589, 828, 589, 828, 590, 829, 590, 830, 590,
830, 591, 831, 591, 831, 592, 832, 592, 832, 593,
833, 593, 834, 593, 834, 594, 835, 594, 835, 595,
836, 595, 837, 595, 837, 596, 838, 596, 838, 597,
839, 597, 840, 597, 840, 598, 841, 598, 841, 599,
842, 599, 842, 600, 843, 600, 844, 600, 844, 601,
845, 601, 845, 602, 846, 602, 847, 602, 847, 603,
848, 603, 848, 604, 849, 604, 849, 605, 850, 605,
851, 605, 851, 606, 852, 606, 852, 607, 853, 607,
854, 607, 854, 608, 855, 608, 855, 609, 856, 609,
856, 610, 857, 610, 858, 610, 858, 611, 859, 611,
859, 612, 860, 612, 861, 612, 861, 613, 862, 613,
862, 614, 863, 614, 863, 615, 864, 615, 865, 615,
865, 616, 866, 616, 866, 617, 867, 617, 868, 617,
868, 618, 869, 618, 869, 619, 870, 619, 871, 619,
871, 620, 872, 620, 872, 621, 873, 621, 873, 622,
874, 622, 875, 622, 875, 623, 876, 623, 876, 624,
877, 624, 878, 624, 878, 625, 879, 625, 879, 626,
880, 626, 880, 627, 881, 627, 882, 627, 882, 628,
883, 628, 883, 629, 884, 629, 885, 629, 885, 630,
886, 630, 886, 631, 887, 631, 887, 632, 888, 632,
889, 632, 889, 633, 890, 633, 890, 634, 891, 634,
892, 634, 892, 635, 893, 635, 893, 636, 894, 636,
894, 637, 895, 637, 896, 637, 896, 638, 897, 638,
897, 639, 898, 639, 899, 639, 899, 640, 900, 640,
900, 641, 901, 641, 902, 641, 902, 642, 903, 642,
903, 643, 904, 643, 904, 644, 905, 644, 906, 644,
906, 645, 907, 645, 907, 646, 908, 646, 909, 646,
909, 647, 910, 647, 910, 648, 911, 648, 911, 649,
912, 649, 913, 649, 913, 650, 914, 650, 914, 651,
915, 651, 916, 651, 916, 652, 917, 652, 917, 653,
918, 653, 918, 654, 919, 654, 920, 654, 920, 655,
921, 655, 921, 656, 922, 656, 923, 656, 923, 657,
924, 657, 924, 658, 925, 658, 925, 659, 926, 659,
927, 659, 927, 660, 928, 660, 928, 661, 929, 661,
930, 661, 930, 662, 931, 662, 931, 663, 932, 663,
933, 663, 933, 664, 934, 664, 934, 665, 935, 665,
935, 666, 936, 666, 937, 666, 937, 667, 938, 667,
938, 668, 939, 668, 940, 668, 940, 669, 941, 669,
941, 670, 942, 670, 942, 671, 943, 671, 944, 671,
944, 672, 945, 672, 945, 673, 946, 673, 947, 673,
947, 674, 948, 674, 948, 675, 949, 675, 949, 676,
950, 676, 951, 676, 951, 677, 952, 677, 952, 678,
953, 678, 954, 678, 954, 679, 955, 679, 955, 680,
956, 680, 956, 681, 957, 681, 958, 681, 958, 682,
959, 682, 959, 683, 960, 683, 961, 683, 961, 684,
962, 684, 962, 685, 963, 685, 964, 685, 964, 686,
965, 686, 965, 687, 966, 687, 966, 688, 967, 688,
968, 688, 968, 689, 969, 689, 969, 690, 970, 690,
971, 690, 971, 691, 972, 691, 972, 692, 973, 692,
973, 693, 974, 693, 975, 693, 975, 694, 976, 694,
976, 695, 977, 695, 978, 695, 978, 696, 979, 696,
979, 697, 980, 697, 980, 698, 981, 698, 982, 698,
982, 699, 983, 699, 983, 700, 984, 700, 985, 700,
985, 701, 986, 701, 986, 702, 987, 702, 988, 702,
988, 703, 989, 703, 989, 704, 990, 704, 990, 705,
991, 705, 992, 705, 992, 706, 993, 706, 993, 707,
994, 707, 995, 707, 995, 708, 996, 708, 996, 709,
997, 709, 997, 710, 998, 710, 999, 710, 999, 711,
1000, 711, 1000, 712, 1001, 712, 1002, 712, 1002, 713,
1003, 713, 1003, 714, 1004, 714, 1004, 715, 1005, 715,
1006, 715, 1006, 716, 1007, 716, 1007, 717, 1008, 717,
1009, 717, 1009, 718, 1010, 718, 1010, 719, 1011, 719,
1011, 720, 1012, 720, 1013, 720, 1013, 721, 1014, 721,
1014, 722, 1015, 722, 1016, 722, 1016, 723, 1017, 723,
1017, 724, 1018, 724, 1019, 724, 1019, 725, 1020, 725,
1020, 726, 1021, 726, 1021, 727, 1022, 727, 1023, 727,
1023, 728, 1024, 728, 1024, 729, 1025, 729, 1026, 729,
1026, 730, 1027, 730, 1027, 731, 1028, 731, 1028, 732,
1029, 732, 1030, 732, 1030, 733, 1031, 733, 1031, 734,
1032, 734, 1033, 734, 1033, 735, 1034, 735, 1034, 736,
1035, 736, 1035, 737, 1036, 737, 1037, 737, 1037, 738,
1038, 738, 1038, 739, 1039, 739, 1040, 739, 1040, 740,
1041, 740, 1041, 741, 1042, 741, 1042, 742, 1043, 742,
1044, 742, 1044, 743, 1045, 743, 1045, 744, 1046, 744,
1047, 744, 1047, 745, 1048, 745, 1048, 746, 1049, 746,
1050, 746, 1050, 747, 1051, 747, 1051, 748, 1052, 748,
1052, 749, 1053, 749, 1054, 749, 1054, 750, 1055, 750,
1055, 751, 1056, 751, 1057, 751, 1057, 752, 1058, 752,
1058, 753, 1059, 753, 1059, 754, 1060, 754, 1061, 754,
1061, 755, 1062, 755, 1062, 756, 1063, 756, 1064, 756,
1064, 757, 1065, 757, 1065, 758, 1066, 758, 1066, 759,
1067, 759, 1068, 759, 1068, 760, 1069, 760, 1069, 761,
1070, 761, 1071, 761, 1071, 762, 1072, 762, 1072, 763,
1073, 763, 1073, 764, 1074, 764, 1075, 764, 1075, 765,
1076, 765, 1076, 766, 1077, 766, 1078, 766, 1078, 767,
1079, 767, 1079, 768, 1080, 768, 1081, 768, 1081, 769,
1082, 769, 1082, 770, 1083, 770, 1083, 771, 1084, 771,
1085, 771, 1085, 772, 1086, 772, 1086, 773, 1087, 773,
1088, 773, 1088, 774, 1089, 774, 1089, 775, 1090, 775,
1090, 776, 1091, 776, 1092, 776, 1092, 777, 1093, 777,
1093, 778, 1094, 778, 1095, 778, 1095, 779, 1096, 779,
1096, 780, 1097, 780, 1097, 781, 1098, 781, 1099, 781,
1099, 782, 1100, 782, 1100, 783, 1101, 783, 1102, 783,
1102, 784, 1103, 784, 1103, 785, 1104, 785, 1104, 786,
1105, 786, 1106, 786, 1106, 787, 1107, 787, 1107, 788,
1108, 788, 1109, 788, 1109, 789, 1110, 789, 1110, 790,
1111, 790, 1112, 790, 1112, 791, 1113, 791, 1113, 792,
1114, 792, 1114, 793, 1115, 793, 1116, 793, 1116, 794,
1117, 794, 1117, 795, 1118, 795, 1119, 795, 1119, 796,
1120, 796, 1120, 797, 1121, 797, 1121, 798, 1122, 798,
1123, 798, 1123, 799, 1124, 799, 1124, 800, 1125, 800,
1126, 800, 1126, 801, 1127, 801, 1127, 802, 1128, 802,
1128, 803, 1129, 803, 1130, 803, 1130, 804, 1131, 804,
1131, 805, 1132, 805, 1133, 805, 1133, 806, 1134, 806,
1134, 807, 1135, 807, 1135, 808, 1136, 808, 1137, 808,
1137, 809, 1138, 809, 1138, 810, 1139, 810, 1140, 810,
1140, 811, 1141, 811, 1141, 812, 1142, 812, 1143, 812,
1143, 813, 1144, 813, 1144, 814, 1145, 814, 1145, 815,
1146, 815, 1147, 815, 1147, 816, 1148, 816, 1148, 817,
1149, 817, 1150, 817, 1150, 818, 1151, 818, 1151, 819,
1152, 819, 1152, 820, 1153, 820, 1154, 820, 1154, 821,
1155, 821, 1155, 822, 1156, 822, 1157, 822, 1157, 823,
1158, 823, 1158, 824, 1159, 824, 1159, 825, 1160, 825,
1161, 825, 1161, 826, 1162, 826, 1162, 827, 1163, 827,
1164, 827, 1164, 828, 1165, 828, 1165, 829, 1166, 829,
1167, 829, 1167, 830, 1168, 830, 1168, 831, 1169, 831,
1169, 832, 1170, 832, 1171, 832, 1171, 833, 1172, 833,
1172, 834, 1173, 834, 1174, 834, 1174, 835, 1175, 835,
1175, 836, 1176, 836, 1176, 837, 1177, 837, 1178, 837,
1178, 838, 1179, 838, 1179, 839, 1180, 839, 1181, 839,
1181, 840, 1182, 840, 1182, 841, 1183, 841, 1183, 842,
1184, 842, 1185, 842, 1185, 843, 1186, 843, 1186, 844,
1187, 844, 1188, 844, 1188, 845, 1189, 845, 1189, 846,
1190, 846, 1190, 847, 1191, 847, 1192, 847, 1192, 848,
1193, 848, 1193, 849, 1194, 849, 1195, 849, 1195, 850,
1196, 850, 1196, 851, 1197, 851, 1198, 851, 1198, 852,
1199, 852, 1199, 853, 1200, 853, 1200, 854, 1201, 854,
1202, 854, 1202, 855, 1203, 855, 1203, 856, 1204, 856,
1205, 856, 1205, 857, 1206, 857, 1206, 858, 1207, 858,
1207, 859, 1208, 859, 1209, 859, 1209, 860, 1210, 860,
1210, 861, 1211, 861, 1212, 861, 1212, 862, 1213, 862,
1213, 863, 1214, 863, 1214, 864, 1215, 864, 1216, 864,
1216, 865, 1217, 865, 1217, 866, 1218, 866, 1219, 866,
1219, 867, 1220, 867, 1220, 868, 1221, 868, 1221, 869,
1222, 869, 1223, 869, 1223, 870, 1224, 870, 1224, 871,
1225, 871, 1226, 871, 1226, 872, 1227, 872, 1227, 873,
1228, 873, 1229, 873, 1229, 874, 1230, 874, 1230, 875,
1231, 875, 1231, 876, 1232, 876, 1233, 876, 1233, 877,
1234, 877, 1234, 878, 1235, 878, 1236, 878, 1236, 879,
1237, 879, 1237, 880, 1238, 880, 1238, 881, 1239, 881,
1240, 881, 1240, 882, 1241, 882, 1241, 883, 1242, 883,
1243, 883, 1243, 884, 1244, 884, 1244, 885, 1245, 885,
1245, 886, 1246, 886, 1247, 886, 1247, 887, 1248, 887,
1248, 888, 1249, 888, 1250, 888, 1250, 889, 1251, 889,
1251, 890, 1252, 890, 1252, 891, 1253, 891, 1254, 891,
1254, 892, 1255, 892, 1255, 893, 1256, 893, 1257, 893,
1257, 894, 1258, 894, 1258, 895, 1259, 895, 1260, 895,
1260, 896, 1261, 896, 1261, 897, 1262, 897, 1262, 898,
1263, 898, 1264, 898, 1264, 899, 1265, 899, 1265, 900,
1266, 900, 1267, 900, 1267, 901, 1268, 901, 1268, 902,
1269, 902, 1269, 903, 1270, 903, 1271, 903, 1271, 904,
1272, 904, 1272, 905, 1273, 905, 1274, 905, 1274, 906,
1275, 906, 1275, 907, 1276, 907, 1276, 908, 1277, 908,
1278, 908, 1278, 909, 1279, 909, 1279, 910, 1280, 910,
1281, 910, 1281, 911, 1282, 911, 1282, 912, 1283, 912,
1283, 913, 1284, 913, 1285, 913, 1285, 914, 1286, 914,
1286, 915, 1287, 915, 1288, 915, 1288, 916, 1289, 916,
1289, 917, 1290, 917, 1291, 917, 1291, 918, 1292, 918,
1292, 919, 1293, 919, 1293, 920, 1294, 920, 1295, 920,
1295, 921, 1296, 921, 1296, 922, 1297, 922, 1298, 922,
1298, 923, 1299, 923, 1299, 924, 1300, 924, 1300, 925,
1301, 925, 1302, 925, 1302, 926, 1303, 926, 1303, 927,
1304, 927, 1305, 927, 1305, 928, 1306, 928, 1306, 929,
1307, 929, 1307, 930, 1308, 930, 1309, 930, 1309, 931,
1310, 931, 1310, 932, 1311, 932, 1312, 932, 1312, 933,
1313, 933, 1313, 934, 1314, 934, 1314, 935, 1315, 935,
1316, 935, 1316, 936, 1317, 936, 1317, 937, 1318, 937,
1319, 937, 1319, 938, 1320, 938, 1320, 939, 1321, 939,
1322, 939, 1322, 940, 1323, 940, 1323, 941, 1324, 941,
1324, 942, 1325, 942, 1326, 942, 1326, 943, 1327, 943,
1327, 944, 1328, 944, 1329, 944, 1329, 945, 1330, 945,
1330, 946, 1331, 946, 1331, 947, 1332, 947, 1333, 947,
1333, 948, 1334, 948, 1334, 949, 1335, 949, 1336, 949,
1336, 950, 1337, 950, 1337, 951, 1338, 951, 1338, 952,
1339, 952, 1340, 952, 1340, 953, 1341, 953, 1341, 954,
1342, 954, 1343, 954, 1343, 955, 1344, 955, 1344, 956,
1345, 956, 1345, 957, 1346, 957, 1347, 957, 1347, 958,
1348, 958, 1348, 959, 1349, 959, 1350, 959, 1350, 960,
1351, 960, 1351, 961, 1352, 961, 1353, 961, 1353, 962,
1354, 962, 1354, 963, 1355, 963, 1355, 964, 1356, 964,
1357, 964, 1357, 965, 1358, 965, 1358, 966, 1359, 966,
1360, 966, 1360, 967, 1361, 967, 1361, 968, 1362, 968,
1362, 969, 1363, 969, 1364, 969, 1364, 970, 1365, 970,
1365, 971, 1366, 971, 1367, 971, 1367, 972, 1368, 972,
1368, 973, 1369, 973, 1369, 974, 1370, 974, 1371, 974,
1371, 975, 1372, 975, 1372, 976, 1373, 976, 1374, 976,
1374, 977, 1375, 977, 1375, 978, 1376, 978, 1377, 978,
1377, 979, 1378, 979, 1378, 980, 1379, 980, 1379, 981,
1380, 981, 1381, 981, 1381, 982, 1382, 982, 1382, 983,
1383, 983, 1384, 983, 1384, 984, 1385, 984, 1385, 985,
1386, 985, 1386, 986, 1387, 986, 1388, 986, 1388, 987,
1389, 987, 1389, 988, 1390, 988, 1391, 988, 1391, 989,
1392, 989, 1392, 990, 1393, 990, 1393, 991, 1394, 991,
1395, 991, 1395, 992, 1396, 992, 1396, 993, 1397, 993,
1398, 993, 1398, 994, 1399, 994, 1399, 995, 1400, 995,
1400, 996, 1401, 996, 1402, 996, 1402, 997, 1403, 997,
1403, 998, 1404, 998, 1405, 998, 1405, 999, 1406, 999,
1406, 1000, 1407, 1000, 1408, 1000, 1408, 1001, 1409, 1001,
1409, 1002, 1410, 1002, 1410, 1003, 1411, 1003, 1412, 1003,
1412, 1004, 1413, 1004, 1413, 1005, 1414, 1005, 1415, 1005,
1415, 1006, 1416, 1006, 1416, 1007, 1417, 1007, 1417, 1008,
1418, 1008, 1419, 1008, 1419, 1009, 1420, 1009, 1420, 1010,
1421, 1010, 1422, 1010, 1422, 1011, 1423, 1011, 1423, 1012,
1424, 1012, 1424, 1013, 1425, 1013, 1426, 1013, 1426, 1014,
1427, 1014, 1427, 1015, 1428, 1015, 1429, 1015, 1429, 1016,
1430, 1016, 1430, 1017, 1431, 1017, 1431, 1018, 1432, 1018,
1433, 1018, 1433, 1019, 1434, 1019, 1434, 1020, 1435, 1020,
1436, 1020, 1436, 1021, 1437, 1021, 1437, 1022, 1438, 1022,
1439, 1022, 1439, 1023, 1440, 1023, 1440, 1024, 1441, 1024,
1441, 1025, 1442, 1025, 1443, 1025, 1443, 1026, 1444, 1026,
1444, 1027, 1445, 1027, 1446, 1027, 1446, 1028, 1447, 1028,
1447, 1029, 1448, 1029, 1448, 1030, 1449, 1030, 1450, 1030,
1450, 1031, 1451, 1031, 1451, 1032, 1452, 1032, 1453, 1032,
1453, 1033, 1454, 1033, 1454, 1034, 1455, 1034, 1455, 1035,
1456, 1035, 1457, 1035, 1457, 1036, 1458, 1036, 1458, 1037,
1459, 1037, 1460, 1037, 1460, 1038, 1461, 1038, 1461, 1039,
1462, 1039, 1462, 1040, 1463, 1040, 1464, 1040, 1464, 1041,
1465, 1041, 1465, 1042, 1466, 1042, 1467, 1042, 1467, 1043,
1468, 1043, 1468, 1044, 1469, 1044, 1470, 1044, 1470, 1045,
1471, 1045, 1471, 1046, 1472, 1046, 1472, 1047, 1473, 1047,
1474, 1047, 1474, 1048, 1475, 1048, 1475, 1049, 1476, 1049,
1477, 1049, 1477, 1050, 1478, 1050, 1478, 1051, 1479, 1051,
1479, 1052, 1480, 1052, 1481, 1052, 1481, 1053, 1482, 1053,
1482, 1054, 1483, 1054, 1484, 1054, 1484, 1055, 1485, 1055,
1485, 1056, 1486, 1056, 1486, 1057, 1487, 1057, 1488, 1057,
1488, 1058, 1489, 1058, 1489, 1059, 1490, 1059, 1491, 1059,
1491, 1060, 1492, 1060, 1492, 1061, 1493, 1061, 1493, 1062,
1494, 1062, 1495, 1062, 1495, 1063, 1496, 1063, 1496, 1064,
1497, 1064, 1498, 1064, 1498, 1065, 1499, 1065, 1499, 1066,
1500, 1066, 1501, 1066, 1501, 1067, 1502, 1067, 1502, 1068,
1503, 1068, 1503, 1069, 1504, 1069, 1505, 1069, 1505, 1070,
1506, 1070, 1506, 1071, 1507, 1071, 1508, 1071, 1508, 1072,
1509, 1072, 1509, 1073};
int main () {
unsigned int m=128; // number of words
int g,h,i,j,k,l,n,a,b,sum;
unsigned int A[256],B[256],C[4],S[256],flag,flag1,length,o;
double x;
FILE *Outfp;
Outfp = fopen("out14b.dat","w");
for (i=0; i<2582; i++) {
l=2*ab[2*i]+ab[2*i+1]+1; // length of loop
n=ab[2*i]+ab[2*i+1]; // number of odd elements in loop
if (n!=(n/d)*d)
continue;
length=3*ab[2*i]+2*ab[2*i+1]; // length of limb
setn(S, 0, m); // clear sum
sum=0; // clear sum
for (j=1; j<=l; j++) {
a=j*n;
if (a==((a/l)*l)) // a=(j*n)/l
a=a/l;
else
a=(a/l)+1;
b=(j-1)*n;
if (b==((b/l)*l)) // b=((j-1)*n)/l
b=b/l;
else
b=(b/l)+1;
if (n<a) {
printf("error \n");
goto zskip;
}
if ((a-b)<0) {
printf("error: a-b is negative \n");
goto zskip;
}
if ((a-b)>0x100000000) {
printf("error: a-b is too large \n");
goto zskip;
}
k=a-b;
k=k<<(j-1); // (a-b)*2**(j-1)
g=1;
setn(A, a-b, m);
for (h=0; h<(n-a); h++) { // (a-b)*3**(n-a)
g=g*3;
copyn(A, B, m);
addn(A, A, m);
addn(B, A, m);
}
lshiftn(A, B, j-1, m); // (a-b)*2**(j-1)*3**(n-a)
addn(B, S, m); // sum+=(a-b)*2**(j-1)*3**(n-a)
if ((S[0]&0xc0000000)!=0) {
printf("error A: sum too large \n");
goto zskip;
}
sum=sum+k*g;
}
if ((sum!=(int)S[m-1])&&(i<18)) {
printf("error: sum=%d \n",sum);
goto zskip;
}
k=1;
k=k<<l;
setn(A, 1, m);
g=1;
for (h=0; h<n; h++) { // 3**n
copyn(A, B, m);
addn(A, A, m);
addn(B, A, m);
if ((A[0]&0xc0000000)!=0) {
printf("error: n=%d \n",n);
goto zskip;
}
g=g*3;
}
if ((unsigned int)l>=(32*m)) {
printf("error: l=%d \n",l);
goto zskip;
}
setn(B, 1, m);
lshiftn(B, B, l, m); // 2**l
subn(B, A, m); // 2**l-3**n
k=k-g;
if ((k!=(int)A[m-1])&&(i<18)) {
printf("error: k=%d \n",k);
goto zskip;
}
x=(double)sum/(double)k;
flag=0;
if ((A[0]&0x80000000)!=0) { // negate 2**l-3**n
setn(B, 0, m);
subn(B, A, m);
flag=1;
}
if ((S[0]&0xff000000)!=0) {
printf("error: sum too big \n");
goto zskip;
}
lshiftn(S, S, 8, m);
flag1=0;
o=orn(A, m-2);
while (o!=0) {
shiftn(S, S, 2, m);
shiftn(A, A, 2, m);
flag1=flag1+2;
o=orn(A, m-2);
}
while ((A[m-2]&0xc0000000)!=0) {
shiftn(S, S, 1, m);
shiftn(A, A, 1, m);
flag1=flag1+1;
}
o=orn(S, m-4);
if ((o|(S[m-4]&0xc0000000))!=0) {
printf("error B: sum too large \n");
goto zskip;
}
div128_64(S[m-4],S[m-3],S[m-2],S[m-1],C,A[m-2],A[m-1]); // sum/(2**l-3**n)
if (i<18)
printf("x=%e \n",x);
fprintf(Outfp,"l=%d, n=%d, length=%d, sign=%d, shift=%d, x=%f \n",l,n,length,flag,flag1,(float)C[3]/256.0);
// fprintf(Outfp,"sum=%#10x %#10x %#10x %#10x \n",S[m-4],S[m-3],S[m-2],S[m-1]);
// fprintf(Outfp,"k=%#10x %#10x \n",A[m-2],A[m-1]);
// fprintf(Outfp,"\n");
shiftn(C, C, 8, 4);
if ((C[0]|C[1]|C[2])!=0) {
printf("warning: x too large \n");
goto zskip;
}
printf("l=%d, n=%d, length=%d, sign=%d, shift=%d, x=%d %d %d %d \n",l,n,length,flag,flag1,C[0],C[1],C[2],C[3]);
}
zskip:
fclose(Outfp);
return(0);
}