n位水仙花数

水仙花数指的是一个 n 位数 ( n≥3 )，它的每个位上的数字的 n 次幂之和等于它本身。(例如:1^3 + 5^3 + 3^3 = 153)

下面代码，求19位，在E7200 2.53G的win xp 32b上耗时1秒左右。欢迎大家指出更多的优化方法。

 

#include <math.h>#include <string.h>#include <stdio.h>#include <stdlib.h>#include <windows.h>unsigned long long g_maxValue;unsigned long long g_minValue;unsigned long long POWXN[10];void preCalculation(unsigned char n){g_maxValue = pow((long double)10, n);g_minValue = pow((long double)10, n-1);POWXN[0] = 0ui64;POWXN[1] = 1ui64;for(int i = 2ui64; i < 10ui64; i++){POWXN[i] = (unsigned long long)(pow((long double)i, n) + 0.5);}}bool search(unsigned long long sum, char *nums, unsigned char n){char *p;int c = 0;int i = 0;int e = n;char tmp;char buf[256];memcpy(buf,nums,n);/*// 19位时，交换花费的计算量仍然超过节约下来的搜索花费。for(i = 0; i < n; ++i){c = sum%10;sum /= 10; p = (char*)memchr(buf,c,e); if(NULL == p) {return false;}else {--e;if(p != buf + e){tmp = buf[e];buf[e] = *p;*p = tmp;}}}*/for(i = 0; i < n; ++i){c = sum%10;sum /= 10; p = (char*)memchr(buf,c,n); if(NULL == p) {return false;}else{*p = 'n';}}return true;}int NarcissusNum(unsigned char n, unsigned long long *out){preCalculation(n);char *nums = (char *)malloc(sizeof(char) * n);memset(nums, 0, sizeof(char) * n);unsigned long long sum = 0;int c = 0;// 遍历所有组合。不用排列，和即是所需排列// 000001, 000002, 000003, ... ,000009, 000011, 000012, ..., 000019, 000022, 000023, ..., 000111, ..., 00122, ..., // 001111, ..., 011111, ..., 111111, ..., 899999, 999999int xx = 0;for(;;){while(nums[0] <= 9){ //xx++;sum = 0;for(unsigned char j = 0; j < n; ++j){sum += POWXN[nums[j]];}if(sum < g_minValue){nums[0]++;continue;}else if(sum > g_maxValue){//break;// 寻找下一个不超过最大值的组合。// 与直接break出去比较，位数越大，优化效果越明显。19位时大约节约10~20msunsigned char i = 1;unsigned char ii;unsigned long long tmp = sum;while(tmp > g_maxValue && i < n){tmp = sum;for(ii = 0; ii < i; ++ii){if(9 == nums[i]){++i;continue;}tmp -= POWXN[nums[ii]];tmp += POWXN[nums[i]+1];}++i;}--i;++nums[i];for(unsigned char ii = 0; ii < i; ++ii){nums[ii] = nums[i];}}else{if(search(sum, nums, n)){*out = sum;out++;}}nums[0]++;}// 组合进位。// 00009 => 000011, 01299 => 01333, 39999 => 44444unsigned char k = 1;for(; k < n; ++k){nums[k]++;if(nums[k] < 10){for(; k > 1; --k){nums[k-1] = nums[k];}break;}}if(nums[n-1] > 9){free(nums);return xx;}nums[0] = nums[1];}free(nums);return xx;} // 校验结果bool Check(unsigned long long val, unsigned char n){unsigned long long tmp = val;unsigned long long sum = 0;unsigned char c;while(0i64 != tmp){c = tmp % 10;tmp /= 10;sum += (unsigned long long )pow((long double)c, n);}if(sum == val)return true;return false;}int _tmain(int argc, _TCHAR* argv[]){ FILETIME ct; FILETIME et; FILETIME kt;FILETIME ut; __int64 ktStart=0; __int64 utStart=0; __int64 ktStop=0;__int64 utStop=0; LARGE_INTEGER ss,ee,fq; HANDLE cur;int i,j,dump;cur = GetCurrentThread(); QueryPerformanceCounter(&ss); //QueryPerformanceCounter的单位应该是cpu的时间周期时间，1/cpu频率GetThreadTimes(cur,&ct,&et,&kt,&ut); // 起始时的内核时间和用户时间。不过精度也就10ms级 ktStart = (((__int64)kt.dwHighDateTime)<<32) + (__int64)kt.dwLowDateTime;utStart = (((__int64)ut.dwHighDateTime)<<32) + (__int64)ut.dwLowDateTime;unsigned long long nums[10] = {0};unsigned char n = 19;for(int c = 0; c < 10; c++){printf("%d/n",NarcissusNum(n, nums));}GetThreadTimes(cur,&ct,&et,&kt,&ut); // 结束时的内核时间和用户时间 ktStop = (((__int64)kt.dwHighDateTime)<<32) + (__int64)kt.dwLowDateTime;utStop = (((__int64)ut.dwHighDateTime)<<32) + (__int64)ut.dwLowDateTime;printf("Kernel Time: %llu ms/n", (ktStop - ktStart)/10000); //FILETIME的单位是100ns，除10000变msprintf("User Time: %llu ms/n", (utStop - utStart)/10000); QueryPerformanceCounter(&ee);QueryPerformanceFrequency(&fq); // QueryPerformanceCounter高精度的，但记得并不是代码实际消耗的时间。printf("Time: %llu us/n", (ee.QuadPart-ss.QuadPart)/(fq.QuadPart/1000000)); for(int i = 0; i < 10 && 0 != nums[i]; i++){if(Check(nums[i], n)){printf("%llu/n", nums[i]);}else{printf("%llu worng/n", nums[i]);}}return 0;}

 

优化了下，800ms

#include <math.h>#include <string.h>#include <stdio.h>#include <stdlib.h>#include <windows.h>unsigned long long g_maxValue;unsigned long long g_minValue;unsigned long long POWXN[10];unsigned long long INC2N[10];void preCalculation(unsigned char n){g_maxValue = pow((long double)10, n);g_minValue = pow((long double)10, n-1);POWXN[0] = 0ui64;POWXN[1] = 1ui64;for(int i = 2; i < 10; i++){POWXN[i] = (unsigned long long)(pow((long double)i, n) + 0.5);}INC2N[0] = 0;INC2N[1] = 1;for(int i = 2; i < 10; i++){INC2N[i] = POWXN[i] - POWXN[i-1];}}bool search(unsigned long long sum, char *nums, unsigned char n){char c;unsigned char nn;unsigned char i;unsigned char j;char buf[256];char *p=buf;memcpy(buf,nums,n);/*for(i = 0; i < n; ++i){c = sum%10;sum /= 10; p = (char*)memchr(buf,c,n); if(NULL == p) {return false;}else{*p = 'n';}}*/if(n < 4){nn = 0;}else{nn = n - 4;}for(i = 0; i < n; i++){c = sum%10;sum /= 10; for(j=0; j < nn; j+=4){if(c == p[j]){p[j] = 'n';c = 'f';break;}else if(c == p[j+1]){p[j+1] = 'n';c = 'f';break;}else if(c == p[j+2]){p[j+2] = 'n';c = 'f';break;}else if(c == p[j+3]){p[j+3] = 'n';c = 'f';break;}} if('f' != c) {for(; j < n; j++){if(c == p[j]){p[j] = 'n';c = 'f';break;}}if('f'!= c)return false;}}return true;}int NarcissusNum(unsigned char n, unsigned long long *out){char c;unsigned char i,j,k,nn;char buf[256];unsigned long long tmp;char *nums;int xx = 0;unsigned long long sum = 0;nums = (char *)malloc(sizeof(char) * n + 1); memset(nums, 0, sizeof(char) * n + 1);memset(buf,0,256);if(n < 4){nn = 0;}else{nn = n - 4;}preCalculation(n);// 遍历所有组合。不用排列，和即是所需排列// 000001, 000002, 000003, ... ,000009, 000011, 000012, ..., 000019, 000022, 000023, ..., 000111, ..., 00122, ..., // 001111, ..., 011111, ..., 111111, ..., 899999, 999999for(;;){while(nums[0] <= 9){ //++xx;sum += INC2N[nums[0]];if(sum < g_minValue){ ++nums[0];continue;}else if(sum > g_maxValue){//++xx; //break;// 寻找下一个不超过最大值的组合。// 与直接break出去比较，位数越大，优化效果越明显。i = 1;tmp = sum;// 假设1 3 8 9 9过大，1 4 4 4 4在范围内while(tmp > g_maxValue && i < n){tmp = sum;// 第一次1 3 9 9 9 仍然过大// 第二次1 4 4 4 4 范围内for(j = 0; j <= i; ++j){if(9 == nums[i]){++i;--j;continue;}tmp -= POWXN[nums[j]];tmp += POWXN[nums[i]+1];}++i;}// tmp 为 1 4 4 4 4的和--i;if(9 == nums[n - 1]){free(nums);return xx;//break;}++nums[i];for(j = 0; j < i; ++j){nums[j] = nums[i];}// sum 降为 1 4 4 4 3的和，上面循环处再加回去。sum = tmp - INC2N[nums[0]];}else{//++xx;/* 19位下，掉用函数大概慢了40msif(search(sum, nums, n)){*out = sum;++out;}++nums[0];*/memcpy(buf,nums,n);tmp = sum;for(i = 0; i < n; i++){c = tmp%10;tmp /= 10; // 展开循环for(j=0; j < nn; j+=4){if(c == buf[j]){buf[j] = 'n';c = 'f';break;}else if(c == buf[j+1]){buf[j+1] = 'n';c = 'f';break;}else if(c == buf[j+2]){buf[j+2] = 'n';c = 'f';break;}else if(c == buf[j+3]){buf[j+3] = 'n';c = 'f';break;}} if('f' != c) {for(; j < n; j++){if(c == buf[j]){buf[j] = 'n';c = 'f';break;}}if('f'!= c){goto NEXT;}}}*out = sum;++out;NEXT:++nums[0];}}// 组合进位。// 00009 => 0000a => 000011, 01299 => 012aa => 01333, 39999 => 3aaaa => 44444//++xx;for(k = 1; k < n; ++k){nums[k]++;if(nums[k] < 10){for(; k > 0; --k){nums[k-1] = nums[k];}break;}}if(nums[n-1] > 9){free(nums);return xx;}sum = POWXN[nums[0]-1]; // 进位后，最低位（nums[0]）不可能是0for(j = 1; j < n; ++j){sum += POWXN[nums[j]];}}free(nums);return xx;} // 校验结果bool Check(unsigned long long val, unsigned char n){unsigned long long tmp = val;unsigned long long sum = 0;unsigned char c;while(0i64 != tmp){c = tmp % 10;tmp /= 10;sum += (unsigned long long )pow((long double)c, n);}if(sum == val)return true;return false;}int _tmain(int argc, _TCHAR* argv[]){ FILETIME ct; FILETIME et; FILETIME kt;FILETIME ut; __int64 ktStart=0; __int64 utStart=0; __int64 ktStop=0;__int64 utStop=0; LARGE_INTEGER ss,ee,fq; HANDLE cur;cur = GetCurrentThread(); QueryPerformanceCounter(&ss); //QueryPerformanceCounter的单位应该是cpu的时间周期时间，1/cpu频率GetThreadTimes(cur,&ct,&et,&kt,&ut); // 起始时的内核时间和用户时间。不过精度也就10ms级 ktStart = (((__int64)kt.dwHighDateTime)<<32) + (__int64)kt.dwLowDateTime;utStart = (((__int64)ut.dwHighDateTime)<<32) + (__int64)ut.dwLowDateTime;unsigned long long nums[10] = {0};unsigned char n = 19;for(int c = 0; c < 10; c++){NarcissusNum(n, nums);//printf("%d/n",NarcissusNum(n, nums));/*NarcissusNum(c, nums);printf("%d:/n",c);for(int i = 0; i < 10 && 0 != nums[i]; i++){if(Check(nums[i], c)){printf("%llu/n", nums[i]);}else{printf("%llu worng/n", nums[i]);}}memset(nums,0,sizeof(unsigned long long) * 10);*/}GetThreadTimes(cur,&ct,&et,&kt,&ut); // 结束时的内核时间和用户时间 ktStop = (((__int64)kt.dwHighDateTime)<<32) + (__int64)kt.dwLowDateTime;utStop = (((__int64)ut.dwHighDateTime)<<32) + (__int64)ut.dwLowDateTime;printf("Kernel Time: %llu ms/n", (ktStop - ktStart)/10000); //FILETIME的单位是100ns，除10000变msprintf("User Time: %llu ms/n", (utStop - utStart)/10000); QueryPerformanceCounter(&ee);QueryPerformanceFrequency(&fq); // QueryPerformanceCounter高精度的，但记得并不是代码实际消耗的时间。printf("Time: %llu us/n", (ee.QuadPart-ss.QuadPart)/(fq.QuadPart/1000000)); //printf("Time: %llu ns/n", (ee.QuadPart-ss.QuadPart)/(fq.QuadPart/1000000000)); for(int i = 0; i < 10 && 0 != nums[i]; i++){if(Check(nums[i], n)){printf("%llu/n", nums[i]);}else{printf("%llu worng/n", nums[i]);}}return 0;}