Cod sursa(job #996719)

Utilizator stefan.petreaStefan Petrea stefan.petrea Data 12 septembrie 2013 15:50:57
Problema Suma divizorilor Scor 40
Compilator cpp Status done
Runda Arhiva de probleme Marime 12.26 kb
#include <iostream>
#include <fstream>
#include <math.h>
using namespace std;
#define min(x,y) (((x)<(y))?(x):(y))
#define max(x,y) (((x)>=(y))?(x):(y))
#define mod 9901
#define pL 1000
// first pL primes
long p[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991, 997, 1009, 1013, 1019, 1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063, 1069, 1087, 1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129, 1151, 1153, 1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217, 1223, 1229, 1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289, 1291, 1297, 1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367, 1373, 1381, 1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447, 1451, 1453, 1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499, 1511, 1523, 1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579, 1583, 1597, 1601, 1607, 1609, 1613, 1619, 1621, 1627, 1637, 1657, 1663, 1667, 1669, 1693, 1697, 1699, 1709, 1721, 1723, 1733, 1741, 1747, 1753, 1759, 1777, 1783, 1787, 1789, 1801, 1811, 1823, 1831, 1847, 1861, 1867, 1871, 1873, 1877, 1879, 1889, 1901, 1907, 1913, 1931, 1933, 1949, 1951, 1973, 1979, 1987, 1993, 1997, 1999, 2003, 2011, 2017, 2027, 2029, 2039, 2053, 2063, 2069, 2081, 2083, 2087, 2089, 2099, 2111, 2113, 2129, 2131, 2137, 2141, 2143, 2153, 2161, 2179, 2203, 2207, 2213, 2221, 2237, 2239, 2243, 2251, 2267, 2269, 2273, 2281, 2287, 2293, 2297, 2309, 2311, 2333, 2339, 2341, 2347, 2351, 2357, 2371, 2377, 2381, 2383, 2389, 2393, 2399, 2411, 2417, 2423, 2437, 2441, 2447, 2459, 2467, 2473, 2477, 2503, 2521, 2531, 2539, 2543, 2549, 2551, 2557, 2579, 2591, 2593, 2609, 2617, 2621, 2633, 2647, 2657, 2659, 2663, 2671, 2677, 2683, 2687, 2689, 2693, 2699, 2707, 2711, 2713, 2719, 2729, 2731, 2741, 2749, 2753, 2767, 2777, 2789, 2791, 2797, 2801, 2803, 2819, 2833, 2837, 2843, 2851, 2857, 2861, 2879, 2887, 2897, 2903, 2909, 2917, 2927, 2939, 2953, 2957, 2963, 2969, 2971, 2999, 3001, 3011, 3019, 3023, 3037, 3041, 3049, 3061, 3067, 3079, 3083, 3089, 3109, 3119, 3121, 3137, 3163, 3167, 3169, 3181, 3187, 3191, 3203, 3209, 3217, 3221, 3229, 3251, 3253, 3257, 3259, 3271, 3299, 3301, 3307, 3313, 3319, 3323, 3329, 3331, 3343, 3347, 3359, 3361, 3371, 3373, 3389, 3391, 3407, 3413, 3433, 3449, 3457, 3461, 3463, 3467, 3469, 3491, 3499, 3511, 3517, 3527, 3529, 3533, 3539, 3541, 3547, 3557, 3559, 3571, 3581, 3583, 3593, 3607, 3613, 3617, 3623, 3631, 3637, 3643, 3659, 3671, 3673, 3677, 3691, 3697, 3701, 3709, 3719, 3727, 3733, 3739, 3761, 3767, 3769, 3779, 3793, 3797, 3803, 3821, 3823, 3833, 3847, 3851, 3853, 3863, 3877, 3881, 3889, 3907, 3911, 3917, 3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989, 4001, 4003, 4007, 4013, 4019, 4021, 4027, 4049, 4051, 4057, 4073, 4079, 4091, 4093, 4099, 4111, 4127, 4129, 4133, 4139, 4153, 4157, 4159, 4177, 4201, 4211, 4217, 4219, 4229, 4231, 4241, 4243, 4253, 4259, 4261, 4271, 4273, 4283, 4289, 4297, 4327, 4337, 4339, 4349, 4357, 4363, 4373, 4391, 4397, 4409, 4421, 4423, 4441, 4447, 4451, 4457, 4463, 4481, 4483, 4493, 4507, 4513, 4517, 4519, 4523, 4547, 4549, 4561, 4567, 4583, 4591, 4597, 4603, 4621, 4637, 4639, 4643, 4649, 4651, 4657, 4663, 4673, 4679, 4691, 4703, 4721, 4723, 4729, 4733, 4751, 4759, 4783, 4787, 4789, 4793, 4799, 4801, 4813, 4817, 4831, 4861, 4871, 4877, 4889, 4903, 4909, 4919, 4931, 4933, 4937, 4943, 4951, 4957, 4967, 4969, 4973, 4987, 4993, 4999, 5003, 5009, 5011, 5021, 5023, 5039, 5051, 5059, 5077, 5081, 5087, 5099, 5101, 5107, 5113, 5119, 5147, 5153, 5167, 5171, 5179, 5189, 5197, 5209, 5227, 5231, 5233, 5237, 5261, 5273, 5279, 5281, 5297, 5303, 5309, 5323, 5333, 5347, 5351, 5381, 5387, 5393, 5399, 5407, 5413, 5417, 5419, 5431, 5437, 5441, 5443, 5449, 5471, 5477, 5479, 5483, 5501, 5503, 5507, 5519, 5521, 5527, 5531, 5557, 5563, 5569, 5573, 5581, 5591, 5623, 5639, 5641, 5647, 5651, 5653, 5657, 5659, 5669, 5683, 5689, 5693, 5701, 5711, 5717, 5737, 5741, 5743, 5749, 5779, 5783, 5791, 5801, 5807, 5813, 5821, 5827, 5839, 5843, 5849, 5851, 5857, 5861, 5867, 5869, 5879, 5881, 5897, 5903, 5923, 5927, 5939, 5953, 5981, 5987, 6007, 6011, 6029, 6037, 6043, 6047, 6053, 6067, 6073, 6079, 6089, 6091, 6101, 6113, 6121, 6131, 6133, 6143, 6151, 6163, 6173, 6197, 6199, 6203, 6211, 6217, 6221, 6229, 6247, 6257, 6263, 6269, 6271, 6277, 6287, 6299, 6301, 6311, 6317, 6323, 6329, 6337, 6343, 6353, 6359, 6361, 6367, 6373, 6379, 6389, 6397, 6421, 6427, 6449, 6451, 6469, 6473, 6481, 6491, 6521, 6529, 6547, 6551, 6553, 6563, 6569, 6571, 6577, 6581, 6599, 6607, 6619, 6637, 6653, 6659, 6661, 6673, 6679, 6689, 6691, 6701, 6703, 6709, 6719, 6733, 6737, 6761, 6763, 6779, 6781, 6791, 6793, 6803, 6823, 6827, 6829, 6833, 6841, 6857, 6863, 6869, 6871, 6883, 6899, 6907, 6911, 6917, 6947, 6949, 6959, 6961, 6967, 6971, 6977, 6983, 6991, 6997, 7001, 7013, 7019, 7027, 7039, 7043, 7057, 7069, 7079, 7103, 7109, 7121, 7127, 7129, 7151, 7159, 7177, 7187, 7193, 7207, 7211, 7213, 7219, 7229, 7237, 7243, 7247, 7253, 7283, 7297, 7307, 7309, 7321, 7331, 7333, 7349, 7351, 7369, 7393, 7411, 7417, 7433, 7451, 7457, 7459, 7477, 7481, 7487, 7489, 7499, 7507, 7517, 7523, 7529, 7537, 7541, 7547, 7549, 7559, 7561, 7573, 7577, 7583, 7589, 7591, 7603, 7607, 7621, 7639, 7643, 7649, 7669, 7673, 7681, 7687, 7691, 7699, 7703, 7717, 7723, 7727, 7741, 7753, 7757, 7759, 7789, 7793, 7817, 7823, 7829, 7841, 7853, 7867, 7873, 7877, 7879, 7883, 7901, 7907, 7919};
long invp[] = {1,4951,7426,8251,8911,9076,6807,9351,9451,8133,9571,9626,4703,5422,5381,8949,6999,9736,9751,1273,4813,5966,5192,4838,6085,9802,2815,6445,6509,9459,8408,7540,9537,5094,5285,9835,2983,1039,7217,2245,5173,9846,469,7993,8234,9851,7025,223,8981,2041,4225,9693,2434,5901,8470,3061,5283,6564,2547,7744,6706,4035,7539,1565,6442,94,9871,3153,9529,6117,6160,4508,2083,9555,6103,959,2118,9876,7401,3179,8314,8463,898,9053,2690,9789,4840,5971,9449,1393,8265,6649,6947,7254,5706,8619,994,1885,4078,3282,6510,5182,2625,889,6757,9265,4545,8747,4222,4934,4395,8720,3520,1522,5642,3388,6832,3050,8428,2874,9886,6527,1450,7012,2999,7058,7202,3213,6287,5992,1677,5164,1967,8002,7543,9424,2347,7621,5759,821,2188,8440,3746,9002,849,7857,3893,5124,6843,5364,1345,2464,4030,5174,7464,3906,2513,8482,8748,3992,8509,9849,4971,8768,4494,8923,9891,2853,1051,4295,8588,5232,8834,8625,9777,8191,1405,5595,1501,9113,3820,4896,3255,7760,3495,4602,3185,1156,6627,7740,9583,1031,3038,2341,4483,3879,2467,7148,5952,4594,2443,8927,761,6561,2369,8869,9196,9399,1253,6900,6318,754,2730,6118,3403,4943,3263,3924,993,6834,6450,2755,1540,1121,4645,6566,1127,9481,8333,676,8094,5616,2418,5479,314,1166,4864,3975,2928,518,618,3488,5965,5750,5480,3921,7960,6237,4436,1283,1706,9276,4334,5831,2205,3078,1094,8945,9435,5375,3068,5628,8328,4418,9685,414,2568,5175,4962,8066,7164,3771,5623,1515,6216,5495,6791,2586,4945,2292,6972,4377,3320,1911,3624,4654,4374,4998,6859,5176,1037,7766,9088,2919,7436,9286,3803,7486,4500,683,6377,1880,2225,7176,6645,5476,3113,7084,2638,2427,581,6315,7897,6672,3631,5970,2034,3280,4327,6230,9478,9507,6010,1077,5498,5771,5396,1857,1333,8500,6953,5504,2650,711,3956,9416,9268,7648,9115,1321,1810,8672,9605,4239,743,1046,1519,7871,3577,3072,4633,1646,7900,5707,5220,6006,7626,6980,2799,7019,7946,5758,7193,8753,4549,5331,6602,5552,6361,1718,9546,7194,1395,4073,4598,265,302,3159,2107,2040,9217,5669,7944,7258,5025,5315,1665,4969,1724,7000,4107,5984,4356,6758,5367,6524,2488,1753,1302,3104,1207,4472,5108,5004,4812,1274,6328,6715,8795,6644,6751,8515,3982,7292,3283,3043,2907,2302,6383,9843,5022,2808,5431,8307,7845,5262,4175,1291,3297,5862,2967,5786,2221,5075,4732,4361,3466,1744,6042,1618,8208,3062,3552,181,8565,9178,5201,3243,6185,1865,7101,7508,1611,4182,7866,8422,365,4262,9423,7891,7874,1632,1414,9898,1165,7638,8967,5636,3698,1335,5703,8022,7669,7246,379,5388,640,3239,7090,798,4186,1341,9795,1943,3103,3048,374,4820,7096,4313,5708,351,5958,4516,8068,7381,3063,142,8095,8567,2431,3480,5058,6867,8122,4893,7613,7198,5683,454,9698,9035,8054,4296,6276,3887,6152,3554,3864,3489,1812,6847,998,7385,5074,9825,1658,1072,5013,1086,930,5444,3133,8484,4439,2274,7109,2192,2106,910,8683,8786,8640,9540,4706,9445,1399,5638,1061,5668,9651,8711,304,4183,3541,8273,6490,3549,4398,3542,2147,522,3774,4469,331,8670,4245,9582,4512,2754,9870,4991,2355,2985,7777,3102,6376,6265,185,9688,3874,7700,8278,4518,4739,8948,9704,9353,949,3676,8704,4592,6078,1783,2698,3348,2978,7057,7178,3349,1969,8421,8427,3405,999,1802,3628,484,634,5323,8765,4132,3532,4258,4679,7490,1931,9508,206,7859,4463,1187,1108,7070,425,3200,2391,2733,6705,5127,2368,871,218,19,7869,8011,9836,9723,3036,2395,8538,2321,8306,3869,1874,6882,3092,2788,1141,8915,3665,9655,7243,6099,2858,3145,9444,8935,7225,1707,5824,9899,8101,1597,5092,8059,4323,9066,1876,8363,6285,9069,1852,2908,6959,6468,3592,4773,6040,8330,8827,6987,3884,2299,532,5343,7739,1654,2526,5723,8296,5989,5633,4079,4428,1020,4491,3311,3326,4429,3711,3969,5873,4542,4664,5502,5132,2404,219,9344,6425,9547,8086,3379,8678,8715,764,1244,981,6275,5254,324,2016,3219,6818,7919,5332,9295,6874,6563,1643,3400,8274,8970,4714,5066,637,6162,8198,3946,7266,3260,7884,6683,5775,4954,3988,1922,7706,7137,3646,4680,2451,9630,5449,3071,329,9560,1286,973,169,9033,1586,8062,3404,965,8696,8511,4568,749,5405,9609,6437,2151,5555,8791,2565,6920,8409,9139,1216,8343,6717,2392,2140,2311,6061,3066,7515,2773,4611,2366,4338,5105,9119,323,3154,872,8155,3021,2871,7726,5590,2456,5815,3008,3568,8034,1776,9008,6580,2224,1714,8410,8522,7815,3690,8548,4438,8485,1442,4394,7331,2796,8043,6724,5853,204,5883,661,8501,8819,6192,5144,8604,359,7977,1380,5050,5720,5224,4026,2823,5768,2113,5667,5493,4353,101,6780,7894,6228,8663,619,7045,9434,2636,6735,5618,2082,2774,9183,5948,6944,2117,5135,8579,3493,4414,145,3539,5453,6570,2667,5054,3885,7888,4350,1234,642,9639,1564,1225,679,2463,4798,5117,962,5644,5725,2265,2461,1487,1791,3478,3905,765,7971,551,3548,2686,682,6113,6805,5924,8024,2979,1024,3381,7485,1854,1862,8299,3849,9138,1285,2002,8373,2372,7682,1262,8384,5159,9538,3496,458,3360,7180,7387,992,4010,5949,8520,2430,4066,2133,4444,8658,6894,2043,1777,8317,4185,2373,8673,8374,6114,1972,6044,8051,5151,4114,2726,2445,2537,130,5743,2045,4885,8721,4502,6392,3501,2624,3873,4789,975,7327,3312,1385,2272,4833,2523,6883,3871,2960,3564,1801,4264,7361,6150,1150,8505,9144,596,3355,7272,7199,5191,4289,1691,9239,3216,9130,1105,6431,8523,2380,6556,7779};

/* the problem actually asks to compute theta_1(A^B) % 9901. 
 *  http://en.wikipedia.org/wiki/Divisor_function
 *
 *
 * theta_1(A^B) =  
 *
 *
 *     r   p_i^{a_i*B} - 1
 * = Prod  ---------------
 *    i=1  p_i - 1
 *
 * because 9901 is prime now, we can compute p_i^{a_i*B} fast through fermat's little theorem
 *
 * we also need the inverse mod 9901 of p_i - 1  which we can also compute fast through extended
 * euclidean algorithm
 *
 */


// powers of primes in A
long a[pL];
int main() {
    ifstream I("sumdiv.in");
    ofstream O("sumdiv.out");
    long A,B;
    int i,j;
    I >> A >> B;

    long aL = -1;//biggest prime divisor index
    for(i=0;i<pL;i++) {
        while(A % p[i] == 0) {
            A /= p[i] ; 
            a[i]++ ; 
            aL = max(i,aL) ; 
            //cout << p[i] << endl;
        };
        if(A == 1)
            break;
    };
    long S = 1;

    long Br = B % (mod-1);
    for(i=0;i<aL+1;i++) {
        if(a[i]) {
            //cout<<"a["<<i<<"]="<<a[i]<<endl;
            //cout<<"p["<<i<<"]="<<p[i]<<endl;
            long u = 1;
            long v = invp[i];
            a[i] %= mod-1;
            long w = ((a[i] * Br) + 1) % mod;
            //cout<<"Br="<<Br<<endl;
            //cout<<"w="<<w<<endl;
            for(j=1;j<=w;j++) {
                u *= p[i];
                u %= mod;
            };
            u = u - 1;
            if(u < 0)
             u += mod;
            //cout<<"u="<<u<<endl;
            //cout<<"v="<<v<<endl;
            S *= (u * v) % mod;
            S %= mod;
        };
    };

    O << S;
    return 0;
};