Greasy Fork

Shined up real nice.

Dozenal OEIS

Gives sequences on OEIS in dozenal, as well as the standard decimal.

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

(function() {
    ////////////////////////////////////////////////////////////////////////////////////////
    // Big Integer Library v. 5.5
    // Created 2000, last modified 2013
    // Leemon Baird
    // www.leemon.com
    //
    // Version history:
    // v 5.5  17 Mar 2013
    //   - two lines of a form like "if (x<0) x+=n" had the "if" changed to "while" to 
    //     handle the case when x<-n. (Thanks to James Ansell for finding that bug)
    // v 5.4  3 Oct 2009
    //   - added "var i" to greaterShift() so i is not global. (Thanks to PŽter Szab— for finding that bug)
    //
    // v 5.3  21 Sep 2009
    //   - added randProbPrime(k) for probable primes
    //   - unrolled loop in mont_ (slightly faster)
    //   - millerRabin now takes a bigInt parameter rather than an int
    //
    // v 5.2  15 Sep 2009
    //   - fixed capitalization in call to int2bigInt in randBigInt
    //     (thanks to Emili Evripidou, Reinhold Behringer, and Samuel Macaleese for finding that bug)
    //
    // v 5.1  8 Oct 2007 
    //   - renamed inverseModInt_ to inverseModInt since it doesn't change its parameters
    //   - added functions GCD and randBigInt, which call GCD_ and randBigInt_
    //   - fixed a bug found by Rob Visser (see comment with his name below)
    //   - improved comments
    //
    // This file is public domain.   You can use it for any purpose without restriction.
    // I do not guarantee that it is correct, so use it at your own risk.  If you use 
    // it for something interesting, I'd appreciate hearing about it.  If you find 
    // any bugs or make any improvements, I'd appreciate hearing about those too.
    // It would also be nice if my name and URL were left in the comments.  But none 
    // of that is required.
    //
    // This code defines a bigInt library for arbitrary-precision integers.
    // A bigInt is an array of integers storing the value in chunks of bpe bits, 
    // little endian (buff[0] is the least significant word).
    // Negative bigInts are stored two's complement.  Almost all the functions treat
    // bigInts as nonnegative.  The few that view them as two's complement say so
    // in their comments.  Some functions assume their parameters have at least one 
    // leading zero element. Functions with an underscore at the end of the name put
    // their answer into one of the arrays passed in, and have unpredictable behavior 
    // in case of overflow, so the caller must make sure the arrays are big enough to 
    // hold the answer.  But the average user should never have to call any of the 
    // underscored functions.  Each important underscored function has a wrapper function 
    // of the same name without the underscore that takes care of the details for you.  
    // For each underscored function where a parameter is modified, that same variable 
    // must not be used as another argument too.  So, you cannot square x by doing 
    // multMod_(x,x,n).  You must use squareMod_(x,n) instead, or do y=dup(x); multMod_(x,y,n).
    // Or simply use the multMod(x,x,n) function without the underscore, where
    // such issues never arise, because non-underscored functions never change
    // their parameters; they always allocate new memory for the answer that is returned.
    //
    // These functions are designed to avoid frequent dynamic memory allocation in the inner loop.
    // For most functions, if it needs a BigInt as a local variable it will actually use
    // a global, and will only allocate to it only when it's not the right size.  This ensures
    // that when a function is called repeatedly with same-sized parameters, it only allocates
    // memory on the first call.
    //
    // Note that for cryptographic purposes, the calls to Math.random() must 
    // be replaced with calls to a better pseudorandom number generator.
    //
    // In the following, "bigInt" means a bigInt with at least one leading zero element,
    // and "integer" means a nonnegative integer less than radix.  In some cases, integer 
    // can be negative.  Negative bigInts are 2s complement.
    // 
    // The following functions do not modify their inputs.
    // Those returning a bigInt, string, or Array will dynamically allocate memory for that value.
    // Those returning a boolean will return the integer 0 (false) or 1 (true).
    // Those returning boolean or int will not allocate memory except possibly on the first 
    // time they're called with a given parameter size.
    // 
    // bigInt  add(x,y)               //return (x+y) for bigInts x and y.  
    // bigInt  addInt(x,n)            //return (x+n) where x is a bigInt and n is an integer.
    // string  bigInt2str(x,base)     //return a string form of bigInt x in a given base, with 2 <= base <= 95
    // int     bitSize(x)             //return how many bits long the bigInt x is, not counting leading zeros
    // bigInt  dup(x)                 //return a copy of bigInt x
    // boolean equals(x,y)            //is the bigInt x equal to the bigint y?
    // boolean equalsInt(x,y)         //is bigint x equal to integer y?
    // bigInt  expand(x,n)            //return a copy of x with at least n elements, adding leading zeros if needed
    // Array   findPrimes(n)          //return array of all primes less than integer n
    // bigInt  GCD(x,y)               //return greatest common divisor of bigInts x and y (each with same number of elements).
    // boolean greater(x,y)           //is x>y?  (x and y are nonnegative bigInts)
    // boolean greaterShift(x,y,shift)//is (x <<(shift*bpe)) > y?
    // bigInt  int2bigInt(t,n,m)      //return a bigInt equal to integer t, with at least n bits and m array elements
    // bigInt  inverseMod(x,n)        //return (x**(-1) mod n) for bigInts x and n.  If no inverse exists, it returns null
    // int     inverseModInt(x,n)     //return x**(-1) mod n, for integers x and n.  Return 0 if there is no inverse
    // boolean isZero(x)              //is the bigInt x equal to zero?
    // boolean millerRabin(x,b)       //does one round of Miller-Rabin base integer b say that bigInt x is possibly prime? (b is bigInt, 1<b<x)
    // boolean millerRabinInt(x,b)    //does one round of Miller-Rabin base integer b say that bigInt x is possibly prime? (b is int,    1<b<x)
    // bigInt  mod(x,n)               //return a new bigInt equal to (x mod n) for bigInts x and n.
    // int     modInt(x,n)            //return x mod n for bigInt x and integer n.
    // bigInt  mult(x,y)              //return x*y for bigInts x and y. This is faster when y<x.
    // bigInt  multMod(x,y,n)         //return (x*y mod n) for bigInts x,y,n.  For greater speed, let y<x.
    // boolean negative(x)            //is bigInt x negative?
    // bigInt  powMod(x,y,n)          //return (x**y mod n) where x,y,n are bigInts and ** is exponentiation.  0**0=1. Faster for odd n.
    // bigInt  randBigInt(n,s)        //return an n-bit random BigInt (n>=1).  If s=1, then the most significant of those n bits is set to 1.
    // bigInt  randTruePrime(k)       //return a new, random, k-bit, true prime bigInt using Maurer's algorithm.
    // bigInt  randProbPrime(k)       //return a new, random, k-bit, probable prime bigInt (probability it's composite less than 2^-80).
    // bigInt  str2bigInt(s,b,n,m)    //return a bigInt for number represented in string s in base b with at least n bits and m array elements
    // bigInt  sub(x,y)               //return (x-y) for bigInts x and y.  Negative answers will be 2s complement
    // bigInt  trim(x,k)              //return a copy of x with exactly k leading zero elements
    //
    //
    // The following functions each have a non-underscored version, which most users should call instead.
    // These functions each write to a single parameter, and the caller is responsible for ensuring the array 
    // passed in is large enough to hold the result. 
    //
    // void    addInt_(x,n)          //do x=x+n where x is a bigInt and n is an integer
    // void    add_(x,y)             //do x=x+y for bigInts x and y
    // void    copy_(x,y)            //do x=y on bigInts x and y
    // void    copyInt_(x,n)         //do x=n on bigInt x and integer n
    // void    GCD_(x,y)             //set x to the greatest common divisor of bigInts x and y, (y is destroyed).  (This never overflows its array).
    // boolean inverseMod_(x,n)      //do x=x**(-1) mod n, for bigInts x and n. Returns 1 (0) if inverse does (doesn't) exist
    // void    mod_(x,n)             //do x=x mod n for bigInts x and n. (This never overflows its array).
    // void    mult_(x,y)            //do x=x*y for bigInts x and y.
    // void    multMod_(x,y,n)       //do x=x*y  mod n for bigInts x,y,n.
    // void    powMod_(x,y,n)        //do x=x**y mod n, where x,y,n are bigInts (n is odd) and ** is exponentiation.  0**0=1.
    // void    randBigInt_(b,n,s)    //do b = an n-bit random BigInt. if s=1, then nth bit (most significant bit) is set to 1. n>=1.
    // void    randTruePrime_(ans,k) //do ans = a random k-bit true random prime (not just probable prime) with 1 in the msb.
    // void    sub_(x,y)             //do x=x-y for bigInts x and y. Negative answers will be 2s complement.
    //
    // The following functions do NOT have a non-underscored version. 
    // They each write a bigInt result to one or more parameters.  The caller is responsible for
    // ensuring the arrays passed in are large enough to hold the results. 
    //
    // void addShift_(x,y,ys)       //do x=x+(y<<(ys*bpe))
    // void carry_(x)               //do carries and borrows so each element of the bigInt x fits in bpe bits.
    // void divide_(x,y,q,r)        //divide x by y giving quotient q and remainder r
    // int  divInt_(x,n)            //do x=floor(x/n) for bigInt x and integer n, and return the remainder. (This never overflows its array).
    // int  eGCD_(x,y,d,a,b)        //sets a,b,d to positive bigInts such that d = GCD_(x,y) = a*x-b*y
    // void halve_(x)               //do x=floor(|x|/2)*sgn(x) for bigInt x in 2's complement.  (This never overflows its array).
    // void leftShift_(x,n)         //left shift bigInt x by n bits.  n<bpe.
    // void linComb_(x,y,a,b)       //do x=a*x+b*y for bigInts x and y and integers a and b
    // void linCombShift_(x,y,b,ys) //do x=x+b*(y<<(ys*bpe)) for bigInts x and y, and integers b and ys
    // void mont_(x,y,n,np)         //Montgomery multiplication (see comments where the function is defined)
    // void multInt_(x,n)           //do x=x*n where x is a bigInt and n is an integer.
    // void rightShift_(x,n)        //right shift bigInt x by n bits.  0 <= n < bpe. (This never overflows its array).
    // void squareMod_(x,n)         //do x=x*x  mod n for bigInts x,n
    // void subShift_(x,y,ys)       //do x=x-(y<<(ys*bpe)). Negative answers will be 2s complement.
    //
    // The following functions are based on algorithms from the _Handbook of Applied Cryptography_
    //    powMod_()           = algorithm 14.94, Montgomery exponentiation
    //    eGCD_,inverseMod_() = algorithm 14.61, Binary extended GCD_
    //    GCD_()              = algorothm 14.57, Lehmer's algorithm
    //    mont_()             = algorithm 14.36, Montgomery multiplication
    //    divide_()           = algorithm 14.20  Multiple-precision division
    //    squareMod_()        = algorithm 14.16  Multiple-precision squaring
    //    randTruePrime_()    = algorithm  4.62, Maurer's algorithm
    //    millerRabin()       = algorithm  4.24, Miller-Rabin algorithm
    //
    // Profiling shows:
    //     randTruePrime_() spends:
    //         10% of its time in calls to powMod_()
    //         85% of its time in calls to millerRabin()
    //     millerRabin() spends:
    //         99% of its time in calls to powMod_()   (always with a base of 2)
    //     powMod_() spends:
    //         94% of its time in calls to mont_()  (almost always with x==y)
    //
    // This suggests there are several ways to speed up this library slightly:
    //     - convert powMod_ to use a Montgomery form of k-ary window (or maybe a Montgomery form of sliding window)
    //         -- this should especially focus on being fast when raising 2 to a power mod n
    //     - convert randTruePrime_() to use a minimum r of 1/3 instead of 1/2 with the appropriate change to the test
    //     - tune the parameters in randTruePrime_(), including c, m, and recLimit
    //     - speed up the single loop in mont_() that takes 95% of the runtime, perhaps by reducing checking
    //       within the loop when all the parameters are the same length.
    //
    // There are several ideas that look like they wouldn't help much at all:
    //     - replacing trial division in randTruePrime_() with a sieve (that speeds up something taking almost no time anyway)
    //     - increase bpe from 15 to 30 (that would help if we had a 32*32->64 multiplier, but not with JavaScript's 32*32->32)
    //     - speeding up mont_(x,y,n,np) when x==y by doing a non-modular, non-Montgomery square
    //       followed by a Montgomery reduction.  The intermediate answer will be twice as long as x, so that
    //       method would be slower.  This is unfortunate because the code currently spends almost all of its time
    //       doing mont_(x,x,...), both for randTruePrime_() and powMod_().  A faster method for Montgomery squaring
    //       would have a large impact on the speed of randTruePrime_() and powMod_().  HAC has a couple of poorly-worded
    //       sentences that seem to imply it's faster to do a non-modular square followed by a single
    //       Montgomery reduction, but that's obviously wrong.
    ////////////////////////////////////////////////////////////////////////////////////////
    
    //globals
    bpe=0;         //bits stored per array element
    mask=0;        //AND this with an array element to chop it down to bpe bits
    radix=mask+1;  //equals 2^bpe.  A single 1 bit to the left of the last bit of mask.
    
    //the digits for converting to different bases
    digitsStr='0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz_=!@#$%^&*()[]{}|;:,.<>/?`~ \\\'\"+-';
    
    //initialize the global variables
    for (bpe=0; (1<<(bpe+1)) > (1<<bpe); bpe++);  //bpe=number of bits in the mantissa on this platform
    bpe>>=1;                   //bpe=number of bits in one element of the array representing the bigInt
    mask=(1<<bpe)-1;           //AND the mask with an integer to get its bpe least significant bits
    radix=mask+1;              //2^bpe.  a single 1 bit to the left of the first bit of mask
    one=int2bigInt(1,1,1);     //constant used in powMod_()
    
    //the following global variables are scratchpad memory to 
    //reduce dynamic memory allocation in the inner loop
    t=new Array(0);
    ss=t;       //used in mult_()
    s0=t;       //used in multMod_(), squareMod_() 
    s1=t;       //used in powMod_(), multMod_(), squareMod_() 
    s2=t;       //used in powMod_(), multMod_()
    s3=t;       //used in powMod_()
    s4=t; s5=t; //used in mod_()
    s6=t;       //used in bigInt2str()
    s7=t;       //used in powMod_()
    T=t;        //used in GCD_()
    sa=t;       //used in mont_()
    mr_x1=t; mr_r=t; mr_a=t;                                      //used in millerRabin()
    eg_v=t; eg_u=t; eg_A=t; eg_B=t; eg_C=t; eg_D=t;               //used in eGCD_(), inverseMod_()
    md_q1=t; md_q2=t; md_q3=t; md_r=t; md_r1=t; md_r2=t; md_tt=t; //used in mod_()
    
    primes=t; pows=t; s_i=t; s_i2=t; s_R=t; s_rm=t; s_q=t; s_n1=t; 
      s_a=t; s_r2=t; s_n=t; s_b=t; s_d=t; s_x1=t; s_x2=t, s_aa=t; //used in randTruePrime_()
      
    rpprb=t; //used in randProbPrimeRounds() (which also uses "primes")
    
    ////////////////////////////////////////////////////////////////////////////////////////
    
    
    //return array of all primes less than integer n
    function findPrimes(n) {
      var i,s,p,ans;
      s=new Array(n);
      for (i=0;i<n;i++)
        s[i]=0;
      s[0]=2;
      p=0;    //first p elements of s are primes, the rest are a sieve
      for(;s[p]<n;) {                  //s[p] is the pth prime
        for(i=s[p]*s[p]; i<n; i+=s[p]) //mark multiples of s[p]
          s[i]=1;
        p++;
        s[p]=s[p-1]+1;
        for(; s[p]<n && s[s[p]]; s[p]++); //find next prime (where s[p]==0)
      }
      ans=new Array(p);
      for(i=0;i<p;i++)
        ans[i]=s[i];
      return ans;
    }
    
    
    //does a single round of Miller-Rabin base b consider x to be a possible prime?
    //x is a bigInt, and b is an integer, with b<x
    function millerRabinInt(x,b) {
      if (mr_x1.length!=x.length) {
        mr_x1=dup(x);
        mr_r=dup(x);
        mr_a=dup(x);
      }
    
      copyInt_(mr_a,b);
      return millerRabin(x,mr_a);
    }
    
    //does a single round of Miller-Rabin base b consider x to be a possible prime?
    //x and b are bigInts with b<x
    function millerRabin(x,b) {
      var i,j,k,s;
    
      if (mr_x1.length!=x.length) {
        mr_x1=dup(x);
        mr_r=dup(x);
        mr_a=dup(x);
      }
    
      copy_(mr_a,b);
      copy_(mr_r,x);
      copy_(mr_x1,x);
    
      addInt_(mr_r,-1);
      addInt_(mr_x1,-1);
    
      //s=the highest power of two that divides mr_r
      k=0;
      for (i=0;i<mr_r.length;i++)
        for (j=1;j<mask;j<<=1)
          if (x[i] & j) {
            s=(k<mr_r.length+bpe ? k : 0); 
             i=mr_r.length;
             j=mask;
          } else
            k++;
    
      if (s)                
        rightShift_(mr_r,s);
    
      powMod_(mr_a,mr_r,x);
    
      if (!equalsInt(mr_a,1) && !equals(mr_a,mr_x1)) {
        j=1;
        while (j<=s-1 && !equals(mr_a,mr_x1)) {
          squareMod_(mr_a,x);
          if (equalsInt(mr_a,1)) {
            return 0;
          }
          j++;
        }
        if (!equals(mr_a,mr_x1)) {
          return 0;
        }
      }
      return 1;  
    }
    
    //returns how many bits long the bigInt is, not counting leading zeros.
    function bitSize(x) {
      var j,z,w;
      for (j=x.length-1; (x[j]==0) && (j>0); j--);
      for (z=0,w=x[j]; w; (w>>=1),z++);
      z+=bpe*j;
      return z;
    }
    
    //return a copy of x with at least n elements, adding leading zeros if needed
    function expand(x,n) {
      var ans=int2bigInt(0,(x.length>n ? x.length : n)*bpe,0);
      copy_(ans,x);
      return ans;
    }
    
    //return a k-bit true random prime using Maurer's algorithm.
    function randTruePrime(k) {
      var ans=int2bigInt(0,k,0);
      randTruePrime_(ans,k);
      return trim(ans,1);
    }
    
    //return a k-bit random probable prime with probability of error < 2^-80
    function randProbPrime(k) {
      if (k>=600) return randProbPrimeRounds(k,2); //numbers from HAC table 4.3
      if (k>=550) return randProbPrimeRounds(k,4);
      if (k>=500) return randProbPrimeRounds(k,5);
      if (k>=400) return randProbPrimeRounds(k,6);
      if (k>=350) return randProbPrimeRounds(k,7);
      if (k>=300) return randProbPrimeRounds(k,9);
      if (k>=250) return randProbPrimeRounds(k,12); //numbers from HAC table 4.4
      if (k>=200) return randProbPrimeRounds(k,15);
      if (k>=150) return randProbPrimeRounds(k,18);
      if (k>=100) return randProbPrimeRounds(k,27);
                  return randProbPrimeRounds(k,40); //number from HAC remark 4.26 (only an estimate)
    }
    
    //return a k-bit probable random prime using n rounds of Miller Rabin (after trial division with small primes)	
    function randProbPrimeRounds(k,n) {
      var ans, i, divisible, B; 
      B=30000;  //B is largest prime to use in trial division
      ans=int2bigInt(0,k,0);
      
      //optimization: try larger and smaller B to find the best limit.
      
      if (primes.length==0)
        primes=findPrimes(30000);  //check for divisibility by primes <=30000
    
      if (rpprb.length!=ans.length)
        rpprb=dup(ans);
    
      for (;;) { //keep trying random values for ans until one appears to be prime
        //optimization: pick a random number times L=2*3*5*...*p, plus a 
        //   random element of the list of all numbers in [0,L) not divisible by any prime up to p.
        //   This can reduce the amount of random number generation.
        
        randBigInt_(ans,k,0); //ans = a random odd number to check
        ans[0] |= 1; 
        divisible=0;
      
        //check ans for divisibility by small primes up to B
        for (i=0; (i<primes.length) && (primes[i]<=B); i++)
          if (modInt(ans,primes[i])==0 && !equalsInt(ans,primes[i])) {
            divisible=1;
            break;
          }      
        
        //optimization: change millerRabin so the base can be bigger than the number being checked, then eliminate the while here.
        
        //do n rounds of Miller Rabin, with random bases less than ans
        for (i=0; i<n && !divisible; i++) {
          randBigInt_(rpprb,k,0);
          while(!greater(ans,rpprb)) //pick a random rpprb that's < ans
            randBigInt_(rpprb,k,0);
          if (!millerRabin(ans,rpprb))
            divisible=1;
        }
        
        if(!divisible)
          return ans;
      }  
    }
    
    //return a new bigInt equal to (x mod n) for bigInts x and n.
    function mod(x,n) {
      var ans=dup(x);
      mod_(ans,n);
      return trim(ans,1);
    }
    
    //return (x+n) where x is a bigInt and n is an integer.
    function addInt(x,n) {
      var ans=expand(x,x.length+1);
      addInt_(ans,n);
      return trim(ans,1);
    }
    
    //return x*y for bigInts x and y. This is faster when y<x.
    function mult(x,y) {
      var ans=expand(x,x.length+y.length);
      mult_(ans,y);
      return trim(ans,1);
    }
    
    //return (x**y mod n) where x,y,n are bigInts and ** is exponentiation.  0**0=1. Faster for odd n.
    function powMod(x,y,n) {
      var ans=expand(x,n.length);  
      powMod_(ans,trim(y,2),trim(n,2),0);  //this should work without the trim, but doesn't
      return trim(ans,1);
    }
    
    //return (x-y) for bigInts x and y.  Negative answers will be 2s complement
    function sub(x,y) {
      var ans=expand(x,(x.length>y.length ? x.length+1 : y.length+1)); 
      sub_(ans,y);
      return trim(ans,1);
    }
    
    //return (x+y) for bigInts x and y.  
    function add(x,y) {
      var ans=expand(x,(x.length>y.length ? x.length+1 : y.length+1)); 
      add_(ans,y);
      return trim(ans,1);
    }
    
    //return (x**(-1) mod n) for bigInts x and n.  If no inverse exists, it returns null
    function inverseMod(x,n) {
      var ans=expand(x,n.length); 
      var s;
      s=inverseMod_(ans,n);
      return s ? trim(ans,1) : null;
    }
    
    //return (x*y mod n) for bigInts x,y,n.  For greater speed, let y<x.
    function multMod(x,y,n) {
      var ans=expand(x,n.length);
      multMod_(ans,y,n);
      return trim(ans,1);
    }
    
    //generate a k-bit true random prime using Maurer's algorithm,
    //and put it into ans.  The bigInt ans must be large enough to hold it.
    function randTruePrime_(ans,k) {
      var c,m,pm,dd,j,r,B,divisible,z,zz,recSize;
    
      if (primes.length==0)
        primes=findPrimes(30000);  //check for divisibility by primes <=30000
    
      if (pows.length==0) {
        pows=new Array(512);
        for (j=0;j<512;j++) {
          pows[j]=Math.pow(2,j/511.-1.);
        }
      }
    
      //c and m should be tuned for a particular machine and value of k, to maximize speed
      c=0.1;  //c=0.1 in HAC
      m=20;   //generate this k-bit number by first recursively generating a number that has between k/2 and k-m bits
      recLimit=20; //stop recursion when k <=recLimit.  Must have recLimit >= 2
    
      if (s_i2.length!=ans.length) {
        s_i2=dup(ans);
        s_R =dup(ans);
        s_n1=dup(ans);
        s_r2=dup(ans);
        s_d =dup(ans);
        s_x1=dup(ans);
        s_x2=dup(ans);
        s_b =dup(ans);
        s_n =dup(ans);
        s_i =dup(ans);
        s_rm=dup(ans);
        s_q =dup(ans);
        s_a =dup(ans);
        s_aa=dup(ans);
      }
    
      if (k <= recLimit) {  //generate small random primes by trial division up to its square root
        pm=(1<<((k+2)>>1))-1; //pm is binary number with all ones, just over sqrt(2^k)
        copyInt_(ans,0);
        for (dd=1;dd;) {
          dd=0;
          ans[0]= 1 | (1<<(k-1)) | Math.floor(Math.random()*(1<<k));  //random, k-bit, odd integer, with msb 1
          for (j=1;(j<primes.length) && ((primes[j]&pm)==primes[j]);j++) { //trial division by all primes 3...sqrt(2^k)
            if (0==(ans[0]%primes[j])) {
              dd=1;
              break;
            }
          }
        }
        carry_(ans);
        return;
      }
    
      B=c*k*k;    //try small primes up to B (or all the primes[] array if the largest is less than B).
      if (k>2*m)  //generate this k-bit number by first recursively generating a number that has between k/2 and k-m bits
        for (r=1; k-k*r<=m; )
          r=pows[Math.floor(Math.random()*512)];   //r=Math.pow(2,Math.random()-1);
      else
        r=.5;
    
      //simulation suggests the more complex algorithm using r=.333 is only slightly faster.
    
      recSize=Math.floor(r*k)+1;
    
      randTruePrime_(s_q,recSize);
      copyInt_(s_i2,0);
      s_i2[Math.floor((k-2)/bpe)] |= (1<<((k-2)%bpe));   //s_i2=2^(k-2)
      divide_(s_i2,s_q,s_i,s_rm);                        //s_i=floor((2^(k-1))/(2q))
    
      z=bitSize(s_i);
    
      for (;;) {
        for (;;) {  //generate z-bit numbers until one falls in the range [0,s_i-1]
          randBigInt_(s_R,z,0);
          if (greater(s_i,s_R))
            break;
        }                //now s_R is in the range [0,s_i-1]
        addInt_(s_R,1);  //now s_R is in the range [1,s_i]
        add_(s_R,s_i);   //now s_R is in the range [s_i+1,2*s_i]
    
        copy_(s_n,s_q);
        mult_(s_n,s_R); 
        multInt_(s_n,2);
        addInt_(s_n,1);    //s_n=2*s_R*s_q+1
        
        copy_(s_r2,s_R);
        multInt_(s_r2,2);  //s_r2=2*s_R
    
        //check s_n for divisibility by small primes up to B
        for (divisible=0,j=0; (j<primes.length) && (primes[j]<B); j++)
          if (modInt(s_n,primes[j])==0 && !equalsInt(s_n,primes[j])) {
            divisible=1;
            break;
          }      
    
        if (!divisible)    //if it passes small primes check, then try a single Miller-Rabin base 2
          if (!millerRabinInt(s_n,2)) //this line represents 75% of the total runtime for randTruePrime_ 
            divisible=1;
    
        if (!divisible) {  //if it passes that test, continue checking s_n
          addInt_(s_n,-3);
          for (j=s_n.length-1;(s_n[j]==0) && (j>0); j--);  //strip leading zeros
          for (zz=0,w=s_n[j]; w; (w>>=1),zz++);
          zz+=bpe*j;                             //zz=number of bits in s_n, ignoring leading zeros
          for (;;) {  //generate z-bit numbers until one falls in the range [0,s_n-1]
            randBigInt_(s_a,zz,0);
            if (greater(s_n,s_a))
              break;
          }                //now s_a is in the range [0,s_n-1]
          addInt_(s_n,3);  //now s_a is in the range [0,s_n-4]
          addInt_(s_a,2);  //now s_a is in the range [2,s_n-2]
          copy_(s_b,s_a);
          copy_(s_n1,s_n);
          addInt_(s_n1,-1);
          powMod_(s_b,s_n1,s_n);   //s_b=s_a^(s_n-1) modulo s_n
          addInt_(s_b,-1);
          if (isZero(s_b)) {
            copy_(s_b,s_a);
            powMod_(s_b,s_r2,s_n);
            addInt_(s_b,-1);
            copy_(s_aa,s_n);
            copy_(s_d,s_b);
            GCD_(s_d,s_n);  //if s_b and s_n are relatively prime, then s_n is a prime
            if (equalsInt(s_d,1)) {
              copy_(ans,s_aa);
              return;     //if we've made it this far, then s_n is absolutely guaranteed to be prime
            }
          }
        }
      }
    }
    
    //Return an n-bit random BigInt (n>=1).  If s=1, then the most significant of those n bits is set to 1.
    function randBigInt(n,s) {
      var a,b;
      a=Math.floor((n-1)/bpe)+2; //# array elements to hold the BigInt with a leading 0 element
      b=int2bigInt(0,0,a);
      randBigInt_(b,n,s);
      return b;
    }
    
    //Set b to an n-bit random BigInt.  If s=1, then the most significant of those n bits is set to 1.
    //Array b must be big enough to hold the result. Must have n>=1
    function randBigInt_(b,n,s) {
      var i,a;
      for (i=0;i<b.length;i++)
        b[i]=0;
      a=Math.floor((n-1)/bpe)+1; //# array elements to hold the BigInt
      for (i=0;i<a;i++) {
        b[i]=Math.floor(Math.random()*(1<<(bpe-1)));
      }
      b[a-1] &= (2<<((n-1)%bpe))-1;
      if (s==1)
        b[a-1] |= (1<<((n-1)%bpe));
    }
    
    //Return the greatest common divisor of bigInts x and y (each with same number of elements).
    function GCD(x,y) {
      var xc,yc;
      xc=dup(x);
      yc=dup(y);
      GCD_(xc,yc);
      return xc;
    }
    
    //set x to the greatest common divisor of bigInts x and y (each with same number of elements).
    //y is destroyed.
    function GCD_(x,y) {
      var i,xp,yp,A,B,C,D,q,sing;
      if (T.length!=x.length)
        T=dup(x);
    
      sing=1;
      while (sing) { //while y has nonzero elements other than y[0]
        sing=0;
        for (i=1;i<y.length;i++) //check if y has nonzero elements other than 0
          if (y[i]) {
            sing=1;
            break;
          }
        if (!sing) break; //quit when y all zero elements except possibly y[0]
    
        for (i=x.length;!x[i] && i>=0;i--);  //find most significant element of x
        xp=x[i];
        yp=y[i];
        A=1; B=0; C=0; D=1;
        while ((yp+C) && (yp+D)) {
          q =Math.floor((xp+A)/(yp+C));
          qp=Math.floor((xp+B)/(yp+D));
          if (q!=qp)
            break;
          t= A-q*C;   A=C;   C=t;    //  do (A,B,xp, C,D,yp) = (C,D,yp, A,B,xp) - q*(0,0,0, C,D,yp)      
          t= B-q*D;   B=D;   D=t;
          t=xp-q*yp; xp=yp; yp=t;
        }
        if (B) {
          copy_(T,x);
          linComb_(x,y,A,B); //x=A*x+B*y
          linComb_(y,T,D,C); //y=D*y+C*T
        } else {
          mod_(x,y);
          copy_(T,x);
          copy_(x,y);
          copy_(y,T);
        } 
      }
      if (y[0]==0)
        return;
      t=modInt(x,y[0]);
      copyInt_(x,y[0]);
      y[0]=t;
      while (y[0]) {
        x[0]%=y[0];
        t=x[0]; x[0]=y[0]; y[0]=t;
      }
    }
    
    //do x=x**(-1) mod n, for bigInts x and n.
    //If no inverse exists, it sets x to zero and returns 0, else it returns 1.
    //The x array must be at least as large as the n array.
    function inverseMod_(x,n) {
      var k=1+2*Math.max(x.length,n.length);
    
      if(!(x[0]&1)  && !(n[0]&1)) {  //if both inputs are even, then inverse doesn't exist
        copyInt_(x,0);
        return 0;
      }
    
      if (eg_u.length!=k) {
        eg_u=new Array(k);
        eg_v=new Array(k);
        eg_A=new Array(k);
        eg_B=new Array(k);
        eg_C=new Array(k);
        eg_D=new Array(k);
      }
    
      copy_(eg_u,x);
      copy_(eg_v,n);
      copyInt_(eg_A,1);
      copyInt_(eg_B,0);
      copyInt_(eg_C,0);
      copyInt_(eg_D,1);
      for (;;) {
        while(!(eg_u[0]&1)) {  //while eg_u is even
          halve_(eg_u);
          if (!(eg_A[0]&1) && !(eg_B[0]&1)) { //if eg_A==eg_B==0 mod 2
            halve_(eg_A);
            halve_(eg_B);      
          } else {
            add_(eg_A,n);  halve_(eg_A);
            sub_(eg_B,x);  halve_(eg_B);
          }
        }
    
        while (!(eg_v[0]&1)) {  //while eg_v is even
          halve_(eg_v);
          if (!(eg_C[0]&1) && !(eg_D[0]&1)) { //if eg_C==eg_D==0 mod 2
            halve_(eg_C);
            halve_(eg_D);      
          } else {
            add_(eg_C,n);  halve_(eg_C);
            sub_(eg_D,x);  halve_(eg_D);
          }
        }
    
        if (!greater(eg_v,eg_u)) { //eg_v <= eg_u
          sub_(eg_u,eg_v);
          sub_(eg_A,eg_C);
          sub_(eg_B,eg_D);
        } else {                   //eg_v > eg_u
          sub_(eg_v,eg_u);
          sub_(eg_C,eg_A);
          sub_(eg_D,eg_B);
        }
      
        if (equalsInt(eg_u,0)) {
          while (negative(eg_C)) //make sure answer is nonnegative
            add_(eg_C,n);
          copy_(x,eg_C);
    
          if (!equalsInt(eg_v,1)) { //if GCD_(x,n)!=1, then there is no inverse
            copyInt_(x,0);
            return 0;
          }
          return 1;
        }
      }
    }
    
    //return x**(-1) mod n, for integers x and n.  Return 0 if there is no inverse
    function inverseModInt(x,n) {
      var a=1,b=0,t;
      for (;;) {
        if (x==1) return a;
        if (x==0) return 0;
        b-=a*Math.floor(n/x);
        n%=x;
    
        if (n==1) return b; //to avoid negatives, change this b to n-b, and each -= to +=
        if (n==0) return 0;
        a-=b*Math.floor(x/n);
        x%=n;
      }
    }
    
    //this deprecated function is for backward compatibility only. 
    function inverseModInt_(x,n) {
       return inverseModInt(x,n);
    }
    
    
    //Given positive bigInts x and y, change the bigints v, a, and b to positive bigInts such that:
    //     v = GCD_(x,y) = a*x-b*y
    //The bigInts v, a, b, must have exactly as many elements as the larger of x and y.
    function eGCD_(x,y,v,a,b) {
      var g=0;
      var k=Math.max(x.length,y.length);
      if (eg_u.length!=k) {
        eg_u=new Array(k);
        eg_A=new Array(k);
        eg_B=new Array(k);
        eg_C=new Array(k);
        eg_D=new Array(k);
      }
      while(!(x[0]&1)  && !(y[0]&1)) {  //while x and y both even
        halve_(x);
        halve_(y);
        g++;
      }
      copy_(eg_u,x);
      copy_(v,y);
      copyInt_(eg_A,1);
      copyInt_(eg_B,0);
      copyInt_(eg_C,0);
      copyInt_(eg_D,1);
      for (;;) {
        while(!(eg_u[0]&1)) {  //while u is even
          halve_(eg_u);
          if (!(eg_A[0]&1) && !(eg_B[0]&1)) { //if A==B==0 mod 2
            halve_(eg_A);
            halve_(eg_B);      
          } else {
            add_(eg_A,y);  halve_(eg_A);
            sub_(eg_B,x);  halve_(eg_B);
          }
        }
    
        while (!(v[0]&1)) {  //while v is even
          halve_(v);
          if (!(eg_C[0]&1) && !(eg_D[0]&1)) { //if C==D==0 mod 2
            halve_(eg_C);
            halve_(eg_D);      
          } else {
            add_(eg_C,y);  halve_(eg_C);
            sub_(eg_D,x);  halve_(eg_D);
          }
        }
    
        if (!greater(v,eg_u)) { //v<=u
          sub_(eg_u,v);
          sub_(eg_A,eg_C);
          sub_(eg_B,eg_D);
        } else {                //v>u
          sub_(v,eg_u);
          sub_(eg_C,eg_A);
          sub_(eg_D,eg_B);
        }
        if (equalsInt(eg_u,0)) {
          while (negative(eg_C)) {   //make sure a (C) is nonnegative
            add_(eg_C,y);
            sub_(eg_D,x);
          }
          multInt_(eg_D,-1);  ///make sure b (D) is nonnegative
          copy_(a,eg_C);
          copy_(b,eg_D);
          leftShift_(v,g);
          return;
        }
      }
    }
    
    
    //is bigInt x negative?
    function negative(x) {
      return ((x[x.length-1]>>(bpe-1))&1);
    }
    
    
    //is (x << (shift*bpe)) > y?
    //x and y are nonnegative bigInts
    //shift is a nonnegative integer
    function greaterShift(x,y,shift) {
      var i, kx=x.length, ky=y.length;
      k=((kx+shift)<ky) ? (kx+shift) : ky;
      for (i=ky-1-shift; i<kx && i>=0; i++) 
        if (x[i]>0)
          return 1; //if there are nonzeros in x to the left of the first column of y, then x is bigger
      for (i=kx-1+shift; i<ky; i++)
        if (y[i]>0)
          return 0; //if there are nonzeros in y to the left of the first column of x, then x is not bigger
      for (i=k-1; i>=shift; i--)
        if      (x[i-shift]>y[i]) return 1;
        else if (x[i-shift]<y[i]) return 0;
      return 0;
    }
    
    //is x > y? (x and y both nonnegative)
    function greater(x,y) {
      var i;
      var k=(x.length<y.length) ? x.length : y.length;
    
      for (i=x.length;i<y.length;i++)
        if (y[i])
          return 0;  //y has more digits
    
      for (i=y.length;i<x.length;i++)
        if (x[i])
          return 1;  //x has more digits
    
      for (i=k-1;i>=0;i--)
        if (x[i]>y[i])
          return 1;
        else if (x[i]<y[i])
          return 0;
      return 0;
    }
    
    //divide x by y giving quotient q and remainder r.  (q=floor(x/y),  r=x mod y).  All 4 are bigints.
    //x must have at least one leading zero element.
    //y must be nonzero.
    //q and r must be arrays that are exactly the same length as x. (Or q can have more).
    //Must have x.length >= y.length >= 2.
    function divide_(x,y,q,r) {
      var kx, ky;
      var i,j,y1,y2,c,a,b;
      copy_(r,x);
      for (ky=y.length;y[ky-1]==0;ky--); //ky is number of elements in y, not including leading zeros
    
      //normalize: ensure the most significant element of y has its highest bit set  
      b=y[ky-1];
      for (a=0; b; a++)
        b>>=1;  
      a=bpe-a;  //a is how many bits to shift so that the high order bit of y is leftmost in its array element
      leftShift_(y,a);  //multiply both by 1<<a now, then divide both by that at the end
      leftShift_(r,a);
    
      //Rob Visser discovered a bug: the following line was originally just before the normalization.
      for (kx=r.length;r[kx-1]==0 && kx>ky;kx--); //kx is number of elements in normalized x, not including leading zeros
    
      copyInt_(q,0);                      // q=0
      while (!greaterShift(y,r,kx-ky)) {  // while (leftShift_(y,kx-ky) <= r) {
        subShift_(r,y,kx-ky);             //   r=r-leftShift_(y,kx-ky)
        q[kx-ky]++;                       //   q[kx-ky]++;
      }                                   // }
    
      for (i=kx-1; i>=ky; i--) {
        if (r[i]==y[ky-1])
          q[i-ky]=mask;
        else
          q[i-ky]=Math.floor((r[i]*radix+r[i-1])/y[ky-1]);	
    
        //The following for(;;) loop is equivalent to the commented while loop, 
        //except that the uncommented version avoids overflow.
        //The commented loop comes from HAC, which assumes r[-1]==y[-1]==0
        //  while (q[i-ky]*(y[ky-1]*radix+y[ky-2]) > r[i]*radix*radix+r[i-1]*radix+r[i-2])
        //    q[i-ky]--;    
        for (;;) {
          y2=(ky>1 ? y[ky-2] : 0)*q[i-ky];
          c=y2>>bpe;
          y2=y2 & mask;
          y1=c+q[i-ky]*y[ky-1];
          c=y1>>bpe;
          y1=y1 & mask;
    
          if (c==r[i] ? y1==r[i-1] ? y2>(i>1 ? r[i-2] : 0) : y1>r[i-1] : c>r[i]) 
            q[i-ky]--;
          else
            break;
        }
    
        linCombShift_(r,y,-q[i-ky],i-ky);    //r=r-q[i-ky]*leftShift_(y,i-ky)
        if (negative(r)) {
          addShift_(r,y,i-ky);         //r=r+leftShift_(y,i-ky)
          q[i-ky]--;
        }
      }
    
      rightShift_(y,a);  //undo the normalization step
      rightShift_(r,a);  //undo the normalization step
    }
    
    //do carries and borrows so each element of the bigInt x fits in bpe bits.
    function carry_(x) {
      var i,k,c,b;
      k=x.length;
      c=0;
      for (i=0;i<k;i++) {
        c+=x[i];
        b=0;
        if (c<0) {
          b=-(c>>bpe);
          c+=b*radix;
        }
        x[i]=c & mask;
        c=(c>>bpe)-b;
      }
    }
    
    //return x mod n for bigInt x and integer n.
    function modInt(x,n) {
      var i,c=0;
      for (i=x.length-1; i>=0; i--)
        c=(c*radix+x[i])%n;
      return c;
    }
    
    //convert the integer t into a bigInt with at least the given number of bits.
    //the returned array stores the bigInt in bpe-bit chunks, little endian (buff[0] is least significant word)
    //Pad the array with leading zeros so that it has at least minSize elements.
    //There will always be at least one leading 0 element.
    function int2bigInt(t,bits,minSize) {   
      var i,k;
      k=Math.ceil(bits/bpe)+1;
      k=minSize>k ? minSize : k;
      buff=new Array(k);
      copyInt_(buff,t);
      return buff;
    }
    
    //return the bigInt given a string representation in a given base.  
    //Pad the array with leading zeros so that it has at least minSize elements.
    //If base=-1, then it reads in a space-separated list of array elements in decimal.
    //The array will always have at least one leading zero, unless base=-1.
    function str2bigInt(s,base,minSize) {
      var d, i, j, x, y, kk;
      var k=s.length;
      if (base==-1) { //comma-separated list of array elements in decimal
        x=new Array(0);
        for (;;) {
          y=new Array(x.length+1);
          for (i=0;i<x.length;i++)
            y[i+1]=x[i];
          y[0]=parseInt(s,10);
          x=y;
          d=s.indexOf(',',0);
          if (d<1) 
            break;
          s=s.substring(d+1);
          if (s.length==0)
            break;
        }
        if (x.length<minSize) {
          y=new Array(minSize);
          copy_(y,x);
          return y;
        }
        return x;
      }
    
      x=int2bigInt(0,base*k,0);
      for (i=0;i<k;i++) {
        d=digitsStr.indexOf(s.substring(i,i+1),0);
        if (base<=36 && d>=36)  //convert lowercase to uppercase if base<=36
          d-=26;
        if (d>=base || d<0) {   //stop at first illegal character
          break;
        }
        multInt_(x,base);
        addInt_(x,d);
      }
    
      for (k=x.length;k>0 && !x[k-1];k--); //strip off leading zeros
      k=minSize>k+1 ? minSize : k+1;
      y=new Array(k);
      kk=k<x.length ? k : x.length;
      for (i=0;i<kk;i++)
        y[i]=x[i];
      for (;i<k;i++)
        y[i]=0;
      return y;
    }
    
    //is bigint x equal to integer y?
    //y must have less than bpe bits
    function equalsInt(x,y) {
      var i;
      if (x[0]!=y)
        return 0;
      for (i=1;i<x.length;i++)
        if (x[i])
          return 0;
      return 1;
    }
    
    //are bigints x and y equal?
    //this works even if x and y are different lengths and have arbitrarily many leading zeros
    function equals(x,y) {
      var i;
      var k=x.length<y.length ? x.length : y.length;
      for (i=0;i<k;i++)
        if (x[i]!=y[i])
          return 0;
      if (x.length>y.length) {
        for (;i<x.length;i++)
          if (x[i])
            return 0;
      } else {
        for (;i<y.length;i++)
          if (y[i])
            return 0;
      }
      return 1;
    }
    
    //is the bigInt x equal to zero?
    function isZero(x) {
      var i;
      for (i=0;i<x.length;i++)
        if (x[i])
          return 0;
      return 1;
    }
    
    //convert a bigInt into a string in a given base, from base 2 up to base 95.
    //Base -1 prints the contents of the array representing the number.
    function bigInt2str(x,base) {
      var i,t,s="";
    
      if (s6.length!=x.length) 
        s6=dup(x);
      else
        copy_(s6,x);
    
      if (base==-1) { //return the list of array contents
        for (i=x.length-1;i>0;i--)
          s+=x[i]+',';
        s+=x[0];
      }
      else { //return it in the given base
        while (!isZero(s6)) {
          t=divInt_(s6,base);  //t=s6 % base; s6=floor(s6/base);
          s=digitsStr.substring(t,t+1)+s;
        }
      }
      if (s.length==0)
        s="0";
      return s;
    }
    
    //returns a duplicate of bigInt x
    function dup(x) {
      var i;
      buff=new Array(x.length);
      copy_(buff,x);
      return buff;
    }
    
    //do x=y on bigInts x and y.  x must be an array at least as big as y (not counting the leading zeros in y).
    function copy_(x,y) {
      var i;
      var k=x.length<y.length ? x.length : y.length;
      for (i=0;i<k;i++)
        x[i]=y[i];
      for (i=k;i<x.length;i++)
        x[i]=0;
    }
    
    //do x=y on bigInt x and integer y.  
    function copyInt_(x,n) {
      var i,c;
      for (c=n,i=0;i<x.length;i++) {
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x+n where x is a bigInt and n is an integer.
    //x must be large enough to hold the result.
    function addInt_(x,n) {
      var i,k,c,b;
      x[0]+=n;
      k=x.length;
      c=0;
      for (i=0;i<k;i++) {
        c+=x[i];
        b=0;
        if (c<0) {
          b=-(c>>bpe);
          c+=b*radix;
        }
        x[i]=c & mask;
        c=(c>>bpe)-b;
        if (!c) return; //stop carrying as soon as the carry is zero
      }
    }
    
    //right shift bigInt x by n bits.  0 <= n < bpe.
    function rightShift_(x,n) {
      var i;
      var k=Math.floor(n/bpe);
      if (k) {
        for (i=0;i<x.length-k;i++) //right shift x by k elements
          x[i]=x[i+k];
        for (;i<x.length;i++)
          x[i]=0;
        n%=bpe;
      }
      for (i=0;i<x.length-1;i++) {
        x[i]=mask & ((x[i+1]<<(bpe-n)) | (x[i]>>n));
      }
      x[i]>>=n;
    }
    
    //do x=floor(|x|/2)*sgn(x) for bigInt x in 2's complement
    function halve_(x) {
      var i;
      for (i=0;i<x.length-1;i++) {
        x[i]=mask & ((x[i+1]<<(bpe-1)) | (x[i]>>1));
      }
      x[i]=(x[i]>>1) | (x[i] & (radix>>1));  //most significant bit stays the same
    }
    
    //left shift bigInt x by n bits.
    function leftShift_(x,n) {
      var i;
      var k=Math.floor(n/bpe);
      if (k) {
        for (i=x.length; i>=k; i--) //left shift x by k elements
          x[i]=x[i-k];
        for (;i>=0;i--)
          x[i]=0;  
        n%=bpe;
      }
      if (!n)
        return;
      for (i=x.length-1;i>0;i--) {
        x[i]=mask & ((x[i]<<n) | (x[i-1]>>(bpe-n)));
      }
      x[i]=mask & (x[i]<<n);
    }
    
    //do x=x*n where x is a bigInt and n is an integer.
    //x must be large enough to hold the result.
    function multInt_(x,n) {
      var i,k,c,b;
      if (!n)
        return;
      k=x.length;
      c=0;
      for (i=0;i<k;i++) {
        c+=x[i]*n;
        b=0;
        if (c<0) {
          b=-(c>>bpe);
          c+=b*radix;
        }
        x[i]=c & mask;
        c=(c>>bpe)-b;
      }
    }
    
    //do x=floor(x/n) for bigInt x and integer n, and return the remainder
    function divInt_(x,n) {
      var i,r=0,s;
      for (i=x.length-1;i>=0;i--) {
        s=r*radix+x[i];
        x[i]=Math.floor(s/n);
        r=s%n;
      }
      return r;
    }
    
    //do the linear combination x=a*x+b*y for bigInts x and y, and integers a and b.
    //x must be large enough to hold the answer.
    function linComb_(x,y,a,b) {
      var i,c,k,kk;
      k=x.length<y.length ? x.length : y.length;
      kk=x.length;
      for (c=0,i=0;i<k;i++) {
        c+=a*x[i]+b*y[i];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;i<kk;i++) {
        c+=a*x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do the linear combination x=a*x+b*(y<<(ys*bpe)) for bigInts x and y, and integers a, b and ys.
    //x must be large enough to hold the answer.
    function linCombShift_(x,y,b,ys) {
      var i,c,k,kk;
      k=x.length<ys+y.length ? x.length : ys+y.length;
      kk=x.length;
      for (c=0,i=ys;i<k;i++) {
        c+=x[i]+b*y[i-ys];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;c && i<kk;i++) {
        c+=x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x+(y<<(ys*bpe)) for bigInts x and y, and integers a,b and ys.
    //x must be large enough to hold the answer.
    function addShift_(x,y,ys) {
      var i,c,k,kk;
      k=x.length<ys+y.length ? x.length : ys+y.length;
      kk=x.length;
      for (c=0,i=ys;i<k;i++) {
        c+=x[i]+y[i-ys];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;c && i<kk;i++) {
        c+=x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x-(y<<(ys*bpe)) for bigInts x and y, and integers a,b and ys.
    //x must be large enough to hold the answer.
    function subShift_(x,y,ys) {
      var i,c,k,kk;
      k=x.length<ys+y.length ? x.length : ys+y.length;
      kk=x.length;
      for (c=0,i=ys;i<k;i++) {
        c+=x[i]-y[i-ys];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;c && i<kk;i++) {
        c+=x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x-y for bigInts x and y.
    //x must be large enough to hold the answer.
    //negative answers will be 2s complement
    function sub_(x,y) {
      var i,c,k,kk;
      k=x.length<y.length ? x.length : y.length;
      for (c=0,i=0;i<k;i++) {
        c+=x[i]-y[i];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;c && i<x.length;i++) {
        c+=x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x+y for bigInts x and y.
    //x must be large enough to hold the answer.
    function add_(x,y) {
      var i,c,k,kk;
      k=x.length<y.length ? x.length : y.length;
      for (c=0,i=0;i<k;i++) {
        c+=x[i]+y[i];
        x[i]=c & mask;
        c>>=bpe;
      }
      for (i=k;c && i<x.length;i++) {
        c+=x[i];
        x[i]=c & mask;
        c>>=bpe;
      }
    }
    
    //do x=x*y for bigInts x and y.  This is faster when y<x.
    function mult_(x,y) {
      var i;
      if (ss.length!=2*x.length)
        ss=new Array(2*x.length);
      copyInt_(ss,0);
      for (i=0;i<y.length;i++)
        if (y[i])
          linCombShift_(ss,x,y[i],i);   //ss=1*ss+y[i]*(x<<(i*bpe))
      copy_(x,ss);
    }
    
    //do x=x mod n for bigInts x and n.
    function mod_(x,n) {
      if (s4.length!=x.length)
        s4=dup(x);
      else
        copy_(s4,x);
      if (s5.length!=x.length)
        s5=dup(x);  
      divide_(s4,n,s5,x);  //x = remainder of s4 / n
    }
    
    //do x=x*y mod n for bigInts x,y,n.
    //for greater speed, let y<x.
    function multMod_(x,y,n) {
      var i;
      if (s0.length!=2*x.length)
        s0=new Array(2*x.length);
      copyInt_(s0,0);
      for (i=0;i<y.length;i++)
        if (y[i])
          linCombShift_(s0,x,y[i],i);   //s0=1*s0+y[i]*(x<<(i*bpe))
      mod_(s0,n);
      copy_(x,s0);
    }
    
    //do x=x*x mod n for bigInts x,n.
    function squareMod_(x,n) {
      var i,j,d,c,kx,kn,k;
      for (kx=x.length; kx>0 && !x[kx-1]; kx--);  //ignore leading zeros in x
      k=kx>n.length ? 2*kx : 2*n.length; //k=# elements in the product, which is twice the elements in the larger of x and n
      if (s0.length!=k) 
        s0=new Array(k);
      copyInt_(s0,0);
      for (i=0;i<kx;i++) {
        c=s0[2*i]+x[i]*x[i];
        s0[2*i]=c & mask;
        c>>=bpe;
        for (j=i+1;j<kx;j++) {
          c=s0[i+j]+2*x[i]*x[j]+c;
          s0[i+j]=(c & mask);
          c>>=bpe;
        }
        s0[i+kx]=c;
      }
      mod_(s0,n);
      copy_(x,s0);
    }
    
    //return x with exactly k leading zero elements
    function trim(x,k) {
      var i,y;
      for (i=x.length; i>0 && !x[i-1]; i--);
      y=new Array(i+k);
      copy_(y,x);
      return y;
    }
    
    //do x=x**y mod n, where x,y,n are bigInts and ** is exponentiation.  0**0=1.
    //this is faster when n is odd.  x usually needs to have as many elements as n.
    function powMod_(x,y,n) {
      var k1,k2,kn,np;
      if(s7.length!=n.length)
        s7=dup(n);
    
      //for even modulus, use a simple square-and-multiply algorithm,
      //rather than using the more complex Montgomery algorithm.
      if ((n[0]&1)==0) {
        copy_(s7,x);
        copyInt_(x,1);
        while(!equalsInt(y,0)) {
          if (y[0]&1)
            multMod_(x,s7,n);
          divInt_(y,2);
          squareMod_(s7,n); 
        }
        return;
      }
    
      //calculate np from n for the Montgomery multiplications
      copyInt_(s7,0);
      for (kn=n.length;kn>0 && !n[kn-1];kn--);
      np=radix-inverseModInt(modInt(n,radix),radix);
      s7[kn]=1;
      multMod_(x ,s7,n);   // x = x * 2**(kn*bp) mod n
    
      if (s3.length!=x.length)
        s3=dup(x);
      else
        copy_(s3,x);
    
      for (k1=y.length-1;k1>0 & !y[k1]; k1--);  //k1=first nonzero element of y
      if (y[k1]==0) {  //anything to the 0th power is 1
        copyInt_(x,1);
        return;
      }
      for (k2=1<<(bpe-1);k2 && !(y[k1] & k2); k2>>=1);  //k2=position of first 1 bit in y[k1]
      for (;;) {
        if (!(k2>>=1)) {  //look at next bit of y
          k1--;
          if (k1<0) {
            mont_(x,one,n,np);
            return;
          }
          k2=1<<(bpe-1);
        }    
        mont_(x,x,n,np);
    
        if (k2 & y[k1]) //if next bit is a 1
          mont_(x,s3,n,np);
      }
    }
    
    
    //do x=x*y*Ri mod n for bigInts x,y,n, 
    //  where Ri = 2**(-kn*bpe) mod n, and kn is the 
    //  number of elements in the n array, not 
    //  counting leading zeros.  
    //x array must have at least as many elemnts as the n array
    //It's OK if x and y are the same variable.
    //must have:
    //  x,y < n
    //  n is odd
    //  np = -(n^(-1)) mod radix
    function mont_(x,y,n,np) {
      var i,j,c,ui,t,ks;
      var kn=n.length;
      var ky=y.length;
    
      if (sa.length!=kn)
        sa=new Array(kn);
        
      copyInt_(sa,0);
    
      for (;kn>0 && n[kn-1]==0;kn--); //ignore leading zeros of n
      for (;ky>0 && y[ky-1]==0;ky--); //ignore leading zeros of y
      ks=sa.length-1; //sa will never have more than this many nonzero elements.  
    
      //the following loop consumes 95% of the runtime for randTruePrime_() and powMod_() for large numbers
      for (i=0; i<kn; i++) {
        t=sa[0]+x[i]*y[0];
        ui=((t & mask) * np) & mask;  //the inner "& mask" was needed on Safari (but not MSIE) at one time
        c=(t+ui*n[0]) >> bpe;
        t=x[i];
        
        //do sa=(sa+x[i]*y+ui*n)/b   where b=2**bpe.  Loop is unrolled 5-fold for speed
        j=1;
        for (;j<ky-4;) { c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++; }    
        for (;j<ky;)   { c+=sa[j]+ui*n[j]+t*y[j];   sa[j-1]=c & mask;   c>>=bpe;   j++; }
        for (;j<kn-4;) { c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++;
                         c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++; }  
        for (;j<kn;)   { c+=sa[j]+ui*n[j];          sa[j-1]=c & mask;   c>>=bpe;   j++; }   
        for (;j<ks;)   { c+=sa[j];                  sa[j-1]=c & mask;   c>>=bpe;   j++; }  
        sa[j-1]=c & mask;
      }
    
      if (!greater(n,sa))
        sub_(sa,n);
      copy_(x,sa);
    }
    
    // BigInt ends
    
    var i;
    var seqElements = Array();
    var els = document.getElementsByTagName('tt');
    for (i = 0; i < els.length; i++)
        if (els[i].parentNode.nodeName == 'TD')
            seqElements.push(els[i]);
    
    for (i = 0; i < seqElements.length; i++) {
        var seqElement = seqElements[i];
        
        var tt = document.createElement('tt');
        tt.innerHTML = convertList('', '', false, seqElement.innerHTML);
        
        var fontElement = seqElement.parentNode.lastElementChild;
        seqElement.parentNode.insertBefore(document.createElement('br'), fontElement);
        seqElement.parentNode.insertBefore(document.createTextNode('Dozenal:'), fontElement);
        seqElement.parentNode.insertBefore(document.createElement('br'), fontElement);
        seqElement.parentNode.insertBefore(tt, fontElement);
        seqElement.parentNode.insertBefore(document.createElement('br'), fontElement);
    }
    
    function convertList(acc, numAcc, inTag, text) {
        if (text == '')
            return (acc + bigInt2str(str2bigInt(numAcc, 10), 12)).replace(/A/g, 'X').replace(/B/g, 'E');
        else {
            var s = text[0], ss = text.substr(1);
            if (inTag) {
                if (s == '>')
                	return convertList(acc + numAcc + '>', '', false, ss);
            }
            else
                if (s == ',' || s == '<' || s == ' ' || s == '-')
                    if (numAcc == '')
                        return convertList(acc + s, '', s == '<', ss);
                    else
                        return convertList(acc + bigInt2str(str2bigInt(numAcc, 10), 12) + s, '', s == '<', ss);
            return convertList(acc, numAcc + s, inTag, ss);
        }
    }
})();