#include "newL2.h" int mynode, totalnodes; /* Algorithm for modular expoentiation.*/ long long ml_exp(long long b, int e, long long m) { long long result = 1; while (e > 0) { /* Multiply in this bits. contribution while using modulus to keep result small*/ /* Note '&' is the bitwise 'and' operator*/ if ((e & 1) == 1) result = (result * b) % m; e >>= 1; b = (b * b) % m; } return result; } /* Another method for modular expoentiation called on integers.*/ int m_expn(int b, int r, int num) { return (int) ml_exp((long long) b, r, (long long) num); } /* Yet another simplified method for expoentiation that may be used when the modulus (n) has already been defined. */ /* The inputs are the base b and the exponent r. */ int m_exp(int b, int r) { return m_expn(b, r, n); } /* Compact method to find total 'Distance to Cycle' and 'Cycle Size' based on previously computed arrays containing results from the individual nodes. */ void computeResults(const int* distToCycle, const int* cycleSize, long long* allToCycleSum, long long* allCycleLengthSum) { long long sumDistToCycle = 0; long long sumCycleSize = 0; int i = 0; for(i = 0; i < n; i++) { sumDistToCycle += distToCycle[i]; sumCycleSize += cycleSize[i]; } *allToCycleSum += sumDistToCycle; *allCycleLengthSum += sumCycleSize; } /* Implementation of the Rabin-Miller test for primality.*/ bool MillerRabin(int num, int k, int q, int a) { int n1 = num-1; if(m_expn(a,q, num) == 1) return true; int i = 0; for(i = 0; i < k; i++) if(m_expn(a,(int)pow(2,i)*q, num) == n1) return true; return false; } /*Test for primiality using an attached array of primes in the header file bn_prime.h.*/ bool isPrime(int num) { int i = 0; for(i = 0; i < 50;i++) { if(primes[i] > (unsigned)num) return true; if((num % primes[i] == 0) && (num!=primes[i])) return false; } int k = 0; int q = num-1; while(q % 2 == 0) { k++; q >>= 1; } srand(time(0)); int a; for(i = 0; i < 10; i++) { a = (rand() % (num-2)) + 1; if(!MillerRabin(num, k, q, a)) return false; } return true; } /* Method to test whether or not an element is a Primitive Root of an inputted modulus. */ bool isPrimRoot(int base) { if(!isPrime(n)) return false; int n_1 = n-1; if ((unsigned)n_1 > (primes[NUMPRIMES-1]*primes[NUMPRIMES-1])) printf("Error in Primitive Root Testing, n could have prime factor too large for testing\n"); int n1 = n_1; int index = 0; int p; while(n1 > 1 && index < NUMPRIMES) { /*find the primes that divide phi(n)*/ if((n1 % primes[index]) == 0) { p = primes[index]; /* divide out that prime all the way so it isn.t tested again*/ while((n1 % primes[index] == 0)) n1/=primes[index]; /*if base^phi(n)/p is 1, not a prim root*/ if(m_exp(base,n_1/p) == 1) return false; /*if(isPrime(n1)) return m_exp(base,n_1/n1) == 1;*/ if(n1 == 50021) return (m_exp(base,n_1/50021) != 1); } index++; } return true; } /*Method to find the greatest common divisors of two elements.*/ int gcd(int a, int b) { if(a== 0) return b; if(b==0) return a; int r = a % b; int d = b; int c; while (r > 0) { c = d; d = r; r = c % d; } return d; } /*Method to determine whether an element is relatively prime to the modulus.*/ bool isRelPrime(int base) { return gcd(base, n-1) == 1; } /*Method to determine the number of integers relatively prime, that is, Euler Phi Function. This is used to find the number of m-ary graphs that will be computers thereby giving the number of trials.*/ double euler(int numb) { double value; int i=0; int result = numb; for(i=2;i*i <= numb;i++) { if (numb % i == 0) result -= result / i; while (numb % i == 0) numb /= i; } if (numb > 1) result -= result / numb; value=result; return value; } /*Method that initializes the arrays based on the size of the modulus.*/ void setArrays(int * cycleSize, bool* visit, int* distToCycle, bool* image) { int i = 0; for(i = 0; i < n; i++) { visit[i] = false; cycleSize[i] = 0; distToCycle[i] = 0; image[i] = false; } } void zeroList(int * listArray) { int i = 0; for(i = 0; i < n; i++) listArray[i] = 0; } void writeTotalResults( int* maxTAll, int* maxCAll, int* terminalAll, int* imageAll, int* allComponents, int* allCyclicNodes, long long* allToCycleSum, long long* allCycleLengthSum, int* oneCycles, int* twoCycles, int* threeCycles, int* fourCycles, int* fiveCycles, int* sevenCycles, int* tenCycles, int* twentyfiveCycles, int* marker, int* orderMarker) /*int ccMinTotal, double expectedComponents, double expectedCyclicNodes, double expectedTerminalNodes, double expectedImageNodes, double expectedTailLength, double expectedCycleLength, double expectedFixedPoints, double expectedMaxCycle, double expectedMaxTail)*/ { double trials=n-3;//euler((n-1)/M_ARY); char fileStr[20]; sprintf(fileStr, "%d_%d_%d.dat", n, M_ARY, mynode); FILE * out = fopen(fileStr, "w"); /* cycles base i sum of cycle size seen from nodes in base i sum of distance to cycle from nodes in base i terminal nodes for base i max cycle for base i max tail for base i cyclic nodes for base i */ double ccMin, numGraphs; int i = 2; for(i = 2; i < n-1; i++) { fprintf(out, "%d %lld %lld %d %d %d %d\n", allComponents[i], allCycleLengthSum[i], allToCycleSum[i], terminalAll[i], maxCAll[i], maxTAll[i], allCyclicNodes[i]); } fclose(out); char fileStr2[20]; sprintf(fileStr2, "distrib_%d_%d.dat", n, M_ARY); FILE * dis = fopen(fileStr2, "w"); for(i=2; i < n-1; i++) { fprintf(dis, "%d %d %d %d %d %d %d %d %d\n", marker[i], oneCycles[i], twoCycles[i], threeCycles[i], fourCycles[i], fiveCycles[i], sevenCycles[i], tenCycles[i], twentyfiveCycles[i]); } fclose(dis); double cComponents = 0; double cComponentsSquared = 0; double cCyclicNodes = 0; double cCyclicNodesSquared = 0; double cImageNodes = 0; /*variance = 0*/ double cMaxCycle = 0; double cMaxCycleSquared = 0; double cMaxTail = 0; double cMaxTailSquared = 0; double cWeightedCycle = 0; double cWeightedCycleSquared = 0; double cWeightedTail = 0; double cWeightedTailSquared = 0; double cFixedPoints = 0; for (i = 2; i < n-1; i++) { cComponents += ((double)allComponents[i]) / trials; cComponentsSquared += (double)allComponents[i] * (double)allComponents[i] / trials; cCyclicNodes += (double)allCyclicNodes[i] / trials; cCyclicNodesSquared += (double)allCyclicNodes[i] * (double)allCyclicNodes[i] / trials; long double cycle = (long double)allCycleLengthSum[i] / (double)(n-1); cWeightedCycle += (long double)(cycle) / trials; cWeightedCycleSquared += (long double)(cycle*cycle) / trials; double tail = (double)allToCycleSum[i] / (double)(n-1); cWeightedTail += tail / trials; cWeightedTailSquared += (double)(tail*tail) / trials; cImageNodes += ((double)(n-1) - terminalAll[i]) / trials; cMaxCycle += (double)maxCAll[i] / trials; cMaxCycleSquared += (double)maxCAll[i]*(double)maxCAll[i] / trials; cMaxTail += (double)maxTAll[i] / trials; cMaxTailSquared += (double)maxTAll[i]*(double)maxTAll[i] / trials; cFixedPoints += (double) oneCycles[i] / trials; } double cComponentsTot = 0; double cComponentsSquaredTot = 0; double cCyclicNodesTot = 0; double cCyclicNodesSquaredTot = 0; double cImageNodesTot = 0; double cMaxCycleTot = 0; double cMaxCycleSquaredTot = 0; double cMaxTailTot = 0; double cMaxTailSquaredTot = 0; double cWeightedCycleTot = 0; double cWeightedCycleSquaredTot = 0; double cWeightedTailTot = 0; double cWeightedTailSquaredTot = 0; MPI_Reduce( &cComponents, &cComponentsTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cComponentsSquared, &cComponentsSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cCyclicNodes, &cCyclicNodesTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cCyclicNodesSquared, &cCyclicNodesSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cImageNodes, &cImageNodesTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cMaxCycle, &cMaxCycleTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cMaxCycleSquared, &cMaxCycleSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cMaxTail, &cMaxTailTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cMaxTailSquared, &cMaxTailSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cWeightedCycle, &cWeightedCycleTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cWeightedCycleSquared, &cWeightedCycleSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cWeightedTail, &cWeightedTailTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); MPI_Reduce( &cWeightedTailSquared, &cWeightedTailSquaredTot, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); if (mynode == 0) { char res[20]; char resAvg[25]; char resExp[25]; char resAll[25]; sprintf(res, "results_%d.dat", n); sprintf(resAvg, "avg_results_%d.csv", n); sprintf(resExp, "exp_results_%d.csv", n); sprintf(resAll, "all_results_%d.csv", n); FILE* r = fopen(res, "w"); FILE* r1 = fopen(resAvg, "w"); FILE* r2 = fopen(resExp, "w"); FILE* r3 = fopen(resAll, "w"); fprintf(r1, "Prime\nOrder, numGraphs, avgCC, avgCyclicN, avgTermN, avgImgN, avgFP, avgTotCycLgth, avgDistToCyc, avgMCyc, avgMTail\n"); fprintf(r2, "Prime\nOrder, numGraphs, expCC, expCyclicN, expTermN, expImgN, expFP, expTotCycLgth, expDistToCyc\n"); fprintf(r3, "Prime\nOrder, numGraphs, avgCC, expCC, avgCyclicN, expCyclicN, avgTermN, expTermN, avgImgN, expImgN, avgFP, expFP, avgTotCycLgth, expTotCycLgth, avgDistToCyc, expDistToCyc, avgMCyc, avgMTail\n"); fprintf(r1, "%d\n", n); fprintf(r2, "%d\n", n); fprintf(r3, "%d\n", n); double ComponentsVariance = cComponentsSquaredTot - cComponentsTot*cComponentsTot; double CyclicNodesVariance = cCyclicNodesSquaredTot - cCyclicNodesTot*cCyclicNodesTot; double WeightedCycleVariance = cWeightedCycleSquaredTot - cWeightedCycleTot*cWeightedCycleTot; double WeightedTailVariance = cWeightedTailSquaredTot - cWeightedTailTot*cWeightedTailTot; double MaxCycleVariance = cMaxCycleSquaredTot - cMaxCycleTot*cMaxCycleTot; double MaxTailVariance = cMaxTailSquaredTot - cMaxTailTot*cMaxTailTot; double avgCC, expCC; double avgCyclicN, expCyclicN; double avgTermN, expTermN; double avgImgN, expImgN; double avgFP, expFP; double avgTotCycLgth, expTotCycLgth; double avgDistToCyc, expDistToCyc; double avgMCyc, expMCyc; double avgMTail, expMTail; for (i=2; i<=(n-1)/2;++i) { if (orderMarker[i] != 0) { ccMin = (n-1)/((double) i); fprintf(r, "Order: %d\n", i); fprintf(r, "Number of graphs: %d\n", orderMarker[i]); fprintf(r, "ccMin: %f\n", ccMin); numGraphs = ((double) orderMarker[i]); fprintf(r1, "%d,%d,", i, orderMarker[i]); fprintf(r2, "%d,%d,", i, orderMarker[i]); fprintf(r3, "%d,%d,", i, orderMarker[i]); avgCC = allComponents[i]/numGraphs; expCC = ccMin*(0.5)*log(i); fprintf(r, "Average connected components: %f\n", avgCC); fprintf(r, "Expected connected components: %f\n", expCC); fprintf(r1, "%f,", avgCC); fprintf(r2, "%f,", expCC); fprintf(r3, "%f,%f,", avgCC, expCC); avgCyclicN = allCyclicNodes[i]/numGraphs; expCyclicN = ccMin*(sqrt((M_PI*i)/2)); fprintf(r, "Average cyclic nodes: %f\n", avgCyclicN); fprintf(r, "Expected cyclic nodes: %f\n", expCyclicN); fprintf(r1, "%f,", avgCyclicN); fprintf(r2, "%f,", expCyclicN); fprintf(r3, "%f,%f,", avgCyclicN, expCyclicN); avgTermN = terminalAll[i]/numGraphs; expTermN = ccMin*((1.0/M_E)*i); fprintf(r, "Average terminal nodes: %f\n", avgTermN); fprintf(r, "Expected terminal nodes: %f\n", expTermN); fprintf(r1, "%f,", avgTermN); fprintf(r2, "%f,", expTermN); fprintf(r3, "%f,%f,", avgTermN, expTermN); avgImgN = imageAll[i]/numGraphs; expImgN = ccMin*((1-(1.0/M_E))*i); fprintf(r, "Average image nodes: %f\n", avgImgN); fprintf(r, "Expected image nodes: %f\n", expImgN); fprintf(r1, "%f,", avgImgN); fprintf(r2, "%f,", expImgN); fprintf(r3, "%f,%f,", avgImgN, expImgN); avgFP = oneCycles[i]/numGraphs; expFP = ccMin*1; fprintf(r, "Average fixed points: %f\n", oneCycles[i]/numGraphs); fprintf(r, "Expected fixed points: %f\n", expFP); fprintf(r1, "%f,", avgFP); fprintf(r2, "%f,", expFP); fprintf(r3, "%f,%f,", avgFP, expFP); avgTotCycLgth = allCycleLengthSum[i]/numGraphs; expTotCycLgth = ccMin*i*(sqrt((M_PI*i)/8)); fprintf(r, "Average total cycle length: %f\n", avgTotCycLgth); fprintf(r, "Expected total cycle length: %f\n", expTotCycLgth); fprintf(r1, "%f,", avgTotCycLgth); fprintf(r2, "%f,", expTotCycLgth); fprintf(r3, "%f,%f,", avgTotCycLgth, expTotCycLgth); avgDistToCyc = allToCycleSum[i]/numGraphs; expDistToCyc = ccMin*i*(sqrt((M_PI*i)/8)); fprintf(r, "Average distance to cycle: %f\n", avgDistToCyc); fprintf(r, "Expected distance to cycle: %f\n", expDistToCyc); fprintf(r1, "%f,", avgDistToCyc); fprintf(r2, "%f\n", expDistToCyc); fprintf(r3, "%f,%f,", avgDistToCyc, expDistToCyc); avgMCyc = maxCAll[i]/numGraphs; avgMTail = maxTAll[i]/numGraphs; fprintf(r, "Average max cycle: %f\n", avgMCyc); fprintf(r, "Average max tail: %f\n\n", avgMTail); fprintf(r1, "%f,%f\n", avgMCyc, avgMTail); fprintf(r3, "%f,%f\n", avgMCyc, avgMTail); } } i = n-1; if (orderMarker[n-1] != 0) { ccMin = (n-1)/((double) i); fprintf(r, "Order: %d\n", i); fprintf(r, "Number of graphs: %d\n", orderMarker[i]); fprintf(r, "ccMin: %f\n", ccMin); numGraphs = ((double) orderMarker[i]); fprintf(r1, "%d,%d,", i, orderMarker[i]); fprintf(r2, "%d,%d,", i, orderMarker[i]); fprintf(r3, "%d,%d,", i, orderMarker[i]); avgCC = allComponents[i]/numGraphs; expCC = ccMin*(0.5)*log(i); fprintf(r, "Average connected components: %f\n", avgCC); fprintf(r, "Expected connected components: %f\n", expCC); fprintf(r1, "%f,", avgCC); fprintf(r2, "%f,", expCC); fprintf(r3, "%f,%f,", avgCC, expCC); avgCyclicN = allCyclicNodes[i]/numGraphs; expCyclicN = ccMin*(sqrt((M_PI*i)/2)); fprintf(r, "Average cyclic nodes: %f\n", avgCyclicN); fprintf(r, "Expected cyclic nodes: %f\n", expCyclicN); fprintf(r1, "%f,", avgCyclicN); fprintf(r2, "%f,", expCyclicN); fprintf(r3, "%f,%f,", avgCyclicN, expCyclicN); avgTermN = terminalAll[i]/numGraphs; expTermN = ccMin*((1.0/M_E)*i); fprintf(r, "Average terminal nodes: %f\n", avgTermN); fprintf(r, "Expected terminal nodes: %f\n", expTermN); fprintf(r1, "%f,", avgTermN); fprintf(r2, "%f,", expTermN); fprintf(r3, "%f,%f,", avgTermN, expTermN); avgImgN = imageAll[i]/numGraphs; expImgN = ccMin*((1-(1.0/M_E))*i); fprintf(r, "Average image nodes: %f\n", avgImgN); fprintf(r, "Expected image nodes: %f\n", expImgN); fprintf(r1, "%f,", avgImgN); fprintf(r2, "%f,", expImgN); fprintf(r3, "%f,%f,", avgImgN, expImgN); avgFP = oneCycles[i]/numGraphs; expFP = ccMin*1; fprintf(r, "Average fixed points: %f\n", avgFP); fprintf(r, "Expected fixed points: %f\n", expFP); fprintf(r1, "%f,", avgFP); fprintf(r2, "%f,", expFP); fprintf(r3, "%f,%f,", avgFP, expFP); avgTotCycLgth = allCycleLengthSum[i]/numGraphs; expTotCycLgth = ccMin*i*(sqrt((M_PI*i)/8)); fprintf(r, "Average total cycle length: %f\n", avgTotCycLgth); fprintf(r, "Expected total cycle length: %f\n", expTotCycLgth); fprintf(r1, "%f,", avgTotCycLgth); fprintf(r2, "%f,", expTotCycLgth); fprintf(r3, "%f,%f,", avgTotCycLgth, expTotCycLgth); avgDistToCyc = allToCycleSum[i]/numGraphs; expDistToCyc = ccMin*i*(sqrt((M_PI*i)/8)); fprintf(r, "Average distance to cycle: %f\n", avgDistToCyc); fprintf(r, "Expected distance to cycle: %f\n", expDistToCyc); fprintf(r1, "%f,", avgDistToCyc); fprintf(r2, "%f\n", expDistToCyc); fprintf(r3, "%f,%f,", avgDistToCyc, expDistToCyc); avgMCyc = maxCAll[i]/numGraphs; avgMTail = maxTAll[i]/numGraphs; fprintf(r, "Average max cycle: %f\n", avgMCyc); fprintf(r, "Average max tail: %f\n\n", avgMTail); fprintf(r1, "%f,%f\n", avgMCyc, avgMTail); fprintf(r3, "%f,%f\n", avgMCyc, avgMTail); } fclose(r); fclose(r1); fclose(r2); fclose(r3); /* //fprintf(r, "components average: %lf \n", cComponentsTot); fprintf(r, "components total: %lf \n", cComponentsTot*trials); fprintf(r, "expected components: %lf \n\n", expectedComponents); //fprintf(r, "components variance: %lf \n", ComponentsVariance); //fprintf(r, "cyclic nodes average: %lf \n", cCyclicNodesTot); fprintf(r, "cyclic nodes total: %lf \n", cCyclicNodesTot*trials); fprintf(r, "expected cyclic nodes: %lf \n\n", expectedCyclicNodes); //fprintf(r, "cyclic nodes variance: %lf\n", CyclicNodesVariance); //fprintf(r, "terminal nodes average: %lf\n",(n-1)-cImageNodesTot); fprintf(r, "terminal nodes total: %lf\n", ((n-1)-cImageNodesTot)*trials); fprintf(r, "expected terminal nodes: %lf\n\n", expectedTerminalNodes); //fprintf(r, "image nodes average: %lf\n", cImageNodesTot); fprintf(r, "image nodes total: %lf\n", cImageNodesTot*trials); fprintf(r, "expected image nodes: %lf\n\n", expectedImageNodes); fprintf(r, "total cycle length (as seen from a node): %lf\n", cWeightedCycle*(trials*(n-1))); fprintf(r, "expected total cycle length: %lf\n\n", expectedCycleLength); fprintf(r, "total distance to a cycle: %lf\n", cWeightedTail*(trials*(n-1))); fprintf(r, "expected total distance to a cycle: %lf\n\n", expectedTailLength); fprintf(r, "avg cycle (as seen from a node): %lf\n", cWeightedCycleTot); fprintf(r, "expected avg cycle: %lf\n\n", expectedCycleLength/ccMinTotal); //fprintf(r, "avg cycle variance: %lf\n", WeightedCycleVariance); fprintf(r, "avg tail (distance from cycle as seen from a node): %lf\n", cWeightedTailTot); fprintf(r, "expected avg tail: %lf\n\n", expectedTailLength/ccMinTotal); //fprintf(r, "avg tail variance: %lf\n", WeightedTailVariance); fprintf(r, "fixed points total: %lf \n",cFixedPoints*trials); fprintf(r, "expected fixed points: %lf \n\n",expectedFixedPoints); fprintf(r, "max cycle average: %lf\n", cMaxCycleTot); fprintf(r, "expected max cycle average: %lf\n\n", expectedMaxCycle); //fprintf(r, "max cycle variance: %lf\n", MaxCycleVariance); fprintf(r, "max tail average: %lf\n", cMaxTailTot); fprintf(r, "expected max tail average: %lf\n\n", expectedMaxTail); //fprintf(r, "max tail variance: %lf\n", MaxTailVariance); */ //fclose(r); //fclose(r1); //fclose(r2); //fclose(r3); } } void run() { MPI_Comm_size(MPI_COMM_WORLD, &totalnodes); MPI_Comm_rank(MPI_COMM_WORLD, &mynode); FILE * s; if (mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "Allocating...\n"); fclose(s); } bool visit[n]; bool image[n]; /*Maximum tail length for base [i]*/ int maxTAll[n]; /*Maximum cycle lenghth for base [i]*/ int maxCAll[n]; /*Terminal nodes for base [i]*/ int terminalAll[n]; int imageAll[n]; /*Size of cycle for the component this node is a part of*/ int cycleSize[n]; /*Distance to cycle from node n (0 if node n is cyclic)*/ int distToCycle[n]; /*Number of components for base i*/ int allComponents[n+1]; /*Number of image nodes for base i*/ int allCyclicNodes[n+1]; /*Sum of all nodes. distance to cycle for base i*/ long long allToCycleSum[n+1]; /*Sum of each node.s cycle length for base i*/ long long allCycleLengthSum[n+1]; /*Record the number of cycles of arbitrary lengths.*/ int oneCycles[n]; int twoCycles[n]; int threeCycles[n]; int fourCycles[n]; int fiveCycles[n]; int sevenCycles[n]; int tenCycles[n]; int twentyfiveCycles[n]; int marker[n]; int orderMarker[n]; /*Initialize variables that will store max cycle length (mC) and max tail length (mT).*/ int mC, mT; int next, loc, baseTail, cycleLength, terminal; int root, exp, base; int listArray[n]; int listSize = 0; int rootOrder, ccMin, ccMinTotal; double expectedComponents=0, expectedCyclicNodes=0, expectedTerminalNodes=0, expectedImageNodes=0, expectedTailLength=0, expectedCycleLength=0; double expectedFixedPoints=0, expectedMaxCycle=0, expectedMaxTail=0; if (mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "zeroing...\n"); fclose(s); } /*initialize arrays to 0 */ int i = 0; for(i = 0; i < n; i++) { maxTAll[i] = 0; maxCAll[i] = 0; terminalAll[i] = 0; imageAll[i] = 0; allComponents[i] = 0; allCyclicNodes[i] = 0; allToCycleSum[i] = 0; allCycleLengthSum[i] = 0; oneCycles[i] = 0; twoCycles[i] = 0; threeCycles[i] = 0; fourCycles[i] = 0; fiveCycles[i] = 0; sevenCycles[i]=0; tenCycles[i] = 0; twentyfiveCycles[i] = 0; marker[i] = 0; orderMarker[i] = 0; } double t; if (mynode == 0) t = MPI_Wtime(); double tt; double expTime = 0; double tailTime = 0; double intoCycleTime = 0; double cycleTime = 0; double resultsTime = 0; if (mynode == 0) { s = fopen(STATUS, "w"); fprintf(s, "Finding a PR...\n"); fclose(s); } //find the smallest primitive root for(root = 1; !isPrimRoot(root); root++); if (mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "Prim root is %d...\n", root); fclose(s); } //int count = -1; for(exp = 1; exp < n-1; exp++) { if (exp == ((n-1)/2)) continue; base = m_exp(root,exp); rootOrder = (n-1)/gcd(exp,n-1); orderMarker[rootOrder] += 1; /* ccMin = (n-1)/rootOrder; ccMinTotal += ccMin; expectedComponents += ccMin*(0.5*log(rootOrder)); expectedCyclicNodes += ccMin*(sqrt((M_PI*rootOrder)/2)); expectedTerminalNodes += ccMin*((1/M_E)*rootOrder); expectedImageNodes += ccMin*((1-(1/M_E))*rootOrder); expectedTailLength += ccMin*(sqrt((M_PI*rootOrder)/8)); expectedCycleLength += ccMin*(sqrt((M_PI*rootOrder)/8)); expectedFixedPoints += ccMin*(1); expectedMaxCycle += ccMin*(0.78248*sqrt(rootOrder)); expectedMaxTail += ccMin*(sqrt(2*M_PI)*log(2)*sqrt(rootOrder));*/ if(exp % 100 == 0 && mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "Exp is %d\n", exp); fclose(s); } //printf("g=%d\n",base); /*discard all but the bases which will make the target M-ARY graphs*/ /*if(gcd(exp,n-1) != M_ARY) continue; count++; if(count % totalnodes != mynode) continue; base = m_exp(root,exp);*/ /*Mark each base that will produce target graph.*/ // marker[base]=1; /*Recall mC is max cycle, mT is max tail.*/ mC = 0; mT = 0; /*0 out everything*/ setArrays(cycleSize, visit, distToCycle, image); /*begin making graph, using gamma(i) = base^i mod n*/ for(i = 1; i < n; i++) { if(visit[i]) continue; next = i; listArray[0] = next; listSize = 1; tt = MPI_Wtime(); while(!visit[next]) { visit[next] = true; next = (next * m_exp(base,next)) % n; image[next] = true; listArray[listSize] = next; listSize++; } expTime += MPI_Wtime() - tt; int j = 0; if(cycleSize[next] != 0) { if(distToCycle[next] == 0) { /*all tail into cycle*/ tt = MPI_Wtime(); cycleLength = cycleSize[listArray[listSize-1]]; if(listSize - 1 > mT) mT = listSize - 1; for(j = 0; j < listSize-1; j++) { distToCycle[listArray[j]] = listSize - 1 - j; cycleSize[listArray[j]] = cycleLength; //printf("All. Node %d, Cycle Length %d\n",listArray[j],cycleLength); } intoCycleTime += MPI_Wtime() - tt; } else { /*extension of tail*/ tt = MPI_Wtime(); baseTail = distToCycle[listArray[listSize-1]]; cycleLength = cycleSize[listArray[listSize-1]]; if(listSize-1 + baseTail > mT) mT = listSize-1 + baseTail; for(j = 0; j < listSize-1; j++) { distToCycle[listArray[j]] = baseTail + listSize - 1 - j; cycleSize[listArray[j]] = cycleLength; //printf("Ext. Node %d, Cycle Length %d\n",listArray[j],cycleLength); } tailTime += MPI_Wtime() - tt; } } else { /*new cycle found*/ tt = MPI_Wtime(); /*loc will be the first node in the cycle we ran in to*/ int repeat = listArray[listSize-1]; for(j = 0; listArray[j] != repeat; j++); int firstCycle = j; cycleLength = listSize - (j+1); if(cycleLength==1) oneCycles[rootOrder]++; if(cycleLength==2) twoCycles[rootOrder]++; if(cycleLength==3) threeCycles[rootOrder]++; if(cycleLength==3) fourCycles[rootOrder]++; if(cycleLength==5) fiveCycles[rootOrder]++; if(cycleLength==7) sevenCycles[rootOrder]++; if(cycleLength==10) tenCycles[rootOrder]++; if(cycleLength==20) twentyfiveCycles[rootOrder]++; if(cycleLength > mC) mC = cycleLength; if(firstCycle > mT) mT = firstCycle; /*mark each tail node along the way with how far it is to the cycle (marked as a negative number)*/ for(j = 0; j < firstCycle; j++) { distToCycle[listArray[j]] = firstCycle - j; cycleSize[listArray[j]] = cycleLength; //printf("Node %d, Cycle Length %d\n",listArray[j],cycleLength); } /*mark each cycle node with how big the cycle is*/ for(j = firstCycle; j < listSize - 1; j++) cycleSize[listArray[j]] = cycleLength; allComponents[rootOrder]++; allCyclicNodes[rootOrder] += cycleLength; cycleTime += MPI_Wtime() - tt; } } tt = MPI_Wtime(); terminal=0; for(i = 1; i < n; i++) if(!image[i]) terminal++; maxTAll[rootOrder] += mT; maxCAll[rootOrder] += mC; terminalAll[rootOrder] += terminal; imageAll[rootOrder] += (n-1)-terminal; computeResults(distToCycle, cycleSize, &allToCycleSum[rootOrder], &allCycleLengthSum[rootOrder]); //printf("all cyle length sum %lld\n",allCycleLengthSum[base]); resultsTime += MPI_Wtime() - tt; } if(mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "Writing Results...\n"); fclose(s); } writeTotalResults( maxTAll, maxCAll, terminalAll, imageAll, allComponents, allCyclicNodes, allToCycleSum, allCycleLengthSum, oneCycles, twoCycles, threeCycles, fourCycles, fiveCycles, sevenCycles, tenCycles, twentyfiveCycles, marker, orderMarker); /*ccMinTotal, expectedComponents, expectedCyclicNodes, expectedTerminalNodes, expectedImageNodes, expectedTailLength, expectedCycleLength, expectedFixedPoints, expectedMaxCycle, expectedMaxTail);*/ if(mynode == 0) { s = fopen(STATUS, "a"); fprintf(s, "%lf minutes...\n Exiting...\n%lf %lf %lf %lf %lf\n", (MPI_Wtime() - t)/60, expTime, tailTime, intoCycleTime, cycleTime, resultsTime); fclose(s); } printf("%lf\n", expTime); }