#include #include #include /* In the following I assume this file has been renamed olson.c Only tested under linux (gcc), but should work with slightly modification on any system. */ /* This is a very simple program for A218030 (for k=2) and similar sequences. Saved as olson.c, it make use of the GMP library for large numbers which must be installed on your system. I compiled it with: gcc -O6 olson.c -lgmp -o olson Invoked as : olson L f2 f3 f5 f7 it finds solutions N < 2^L to : k*N = product of nonzero digits of N^2, where k = 2^f2 * 3^f3 * 5^f5 * 7^f7 and the f may be negative. For example, for k = 35/8, with N<2^100, the command is: olsen 100 -3 0 1 1 The values found are saved in a file named Ols-k-L.txt, where k can be of the form x:y. Note the values of N are usually not in order, so they must be sorted using some other program. NOTE: It could be optimized. It is not to be taken as an example of good programming, because it is not. Since I'm lazy, I did not implement any check on the consistency of the input. Giovanni Resta, Oct 20 2012. */ void Usage() { printf("\n---------------------------------------\n"); printf("USAGE: olson L f2 f3 f5 f7\n\n"); printf("it finds solutions N < 2^L to :\n"); printf("k*N = product of nonzero digits of N^2,\n"); printf("where k = 2^f2 * 3^f3 * 5^f5 * 7^f7\n"); printf("Example, for k = 35/8 and N<2^100, it is\n"); printf("olson 100 -3 0 1 1\n\n"); printf("N are saved in a file names Ols-k-L.txt\n"); printf("and are not sorted\n"); printf("---------------------------------------\n\n"); exit(0); } // used to create the file name void setK(int v, int m, int *num, int *den) { int i; if (v>0) for (i=1; i<=v; i++) *num*=m; else for (i=-1; i>=v; i--) *den*=m; } int main(int argc, char *argv[]) { // the search is limited to numbers below 2^Lim unsigned Lim; // exponents in the factorization of current N int e2, e3, e5, e7; // exponents in the product of the digit of N^2 int c2, c3, c5, c7; // exponents of 2,3,5 and 7 in the factor K // in the solution of // K*N = product of nonzero digits of N^2 int f2, f3, f5, f7; // z2,z3,z5,zN will contain partial products // zLim = 2^Lim is the search limit // zQ holds the square of N mpz_t z2,z3,z5,zN, zLim, zQ; // it holds the actual digits of N^2 char *S; // name of the output file and output file char fname[64]; FILE *f; // working vars for generating the output file name int num,den; // solutions found until now int nsol=0; // ================== START ========== // not the input I expected, print usage and exit if (argc!=6) Usage(); // input the input Lim = atoi(argv[1]); f2 = atoi(argv[2]); f3 = atoi(argv[3]); f5 = atoi(argv[4]); f7 = atoi(argv[5]); // TODO: (not by me... )sanity check on the input values // I build the name of the file num=den=1; setK(f2,2,&num,&den); setK(f3,3,&num,&den); setK(f5,5,&num,&den); setK(f7,7,&num,&den); if (den==1) sprintf(fname,"Ols-%d-%d.txt",num,Lim); else sprintf(fname,"Ols-%d:%d-%d.txt",num,den,Lim); // I reserve a space to write the square digits S = malloc(2+(size_t)(0.301029995664*Lim*2)); // init the GMP variables mpz_init(z2); mpz_init(z3); mpz_init(z5); mpz_init(zLim); mpz_init(zN); mpz_init(zQ); // set the search limit as zLim=2^Lim mpz_ui_pow_ui(zLim,2,Lim); // search all the numbers 2^e2*3^e3*5^e5*7^e7 < Lim mpz_set_ui(z2,1); e2=0; while (e2 2^%u [found = %d]\n",e2,nsol); mpz_set(z3,z2); e3=0; while (mpz_cmp(z3,zLim)<0) { mpz_set(z5,z3); e5=0; while (mpz_cmp(z5,zLim)<0) { mpz_set(zN,z5); e7=0; while (mpz_cmp(zN,zLim)<0) { char *p; c2=c3=c5=c7=0; // Inner loop, now I have my N and I check it //-------------------- // zQ = zN^2 mpz_mul(zQ,zN,zN); // S = string corresponding to zQ=zN^2 mpz_get_str(S,10,zQ); /* I scan the digits of N^2 and I update the exponents of the factorization of their product */ for (p=S; *p; p++) switch (*p) { case '2': c2++; break; case '3': c3++; break; case '4': c2+=2; break; case '5': c5++; break; case '6': c2++; c3++; break; case '7': c7++; break; case '8': c2+=3; break; case '9': c3+=2; break; } // I check if, based of their factors // the product of digits is k*N if (c2==f2+e2 && c3==f3+e3 && c5==f5+e5 && c7==f7+e7) { // I output the solution found nsol++; mpz_get_str(S,10,zN); printf("%s\n",S); f=fopen(fname,"a"); fprintf(f,"%s\n",S); fclose(f); } //-------------------- mpz_mul_ui(zN,zN,7); e7++; } mpz_mul_ui(z5,z5,5); e5++; } mpz_mul_ui(z3,z3,3); e3++; } mpz_mul_ui(z2,z2,2); e2++; } // I add a line to output file to certify I // ended the search f=fopen(fname,"a"); fprintf(f,"# Search for n < 2^%d completed\n",Lim); fclose(f); return 0; }