求形如 2^r*3^s*5^t 的整数序列

mathematica · 发表于 2013-10-22 09:21:44

gxqcn · 发表于 2013-10-22 09:39:16

mathematica 发表于 2013-10-22 09:23
讨论的是别人剩下来的东西

请注意你的措词！
什么是“别人剩下来的东西”？！

首先，问题本身是个纯数学问题，不可能是什么新问题；但作为算法题，大家来PK还是很不错的（这才是关键）！
其次，发这个问题之前，我并不知道别处是否有过类似讨论，因为我是将数年前在CSDN上看到的题目再推广的。

发这个帖子的目的是为了追求极致的算法，即如何快速解决问题的过程，而非解决问题本身。
所以，何来“别人剩下来的东西”之说？

mathematica · 发表于 2013-10-22 09:44:18

算法求解: 2^i * 3^j * 5^k i,j,k >=0. 如何用程序打印出下列的数 1 2 3 4 5 6 8 9 10 12 15 16 18 20 ..

http://zhidao.baidu.com/question/299302012

mathematica · 发表于 2013-10-22 14:54:17

http://acm.fjnu.edu.cn/problem/1426
#include<stdio.h>
#include<iostream.h>
int min(int a,int b)
{
return (a>b?b:a);
}
int getmin(int a,int b,int c)
{
return min(min(a,b),c);
}
void f()
{
int i2,i3,i5,i,a2,a3,a5,n;
while(scanf("%d%d%d%d",&a2,&a3,&a5,&n)!=EOF)
{
int *a=new int[n+1];
a[0]=1;
i2=i5=i3=0;
for(i=1;i<=n;i++)
{
a[ i ]=getmin(a2*a[i2],a3*a[i3],a5*a[i5]);
if(a[ i ]==a2*a[i2])
i2++;
if(a[ i ]==a3*a[i3])
i3++;
if(a[ i ]==a5*a[i5])
i5++;
}
printf("%d\n",a[n]);
delete a;
}
}
int main()
{
f();
return 1;
}

复制代码

我把代码复制过来,相当于备份一下
http://zhidao.baidu.com/question/299302012

wayne · 发表于 2013-10-22 21:25:16

northwolves 发表于 2013-10-21 22:21
这个链接http://rosettacode.org/wiki/Hamming_numbers#PARI.2FGP 的代码：

(22:18) gp >

神链接。

#include <iostream>
#include <vector>
// Hamming like sequences Generator
//
// Nigel Galloway. August 13th., 2012
//
class Ham {
private:
std::vector<unsigned int> _H, _hp, _hv, _x;
public:
bool operator!=(const Ham& other) const {return true;}
Ham begin() const {return *this;}
Ham end() const {return *this;}
unsigned int operator*() const {return _x.back();}
Ham(const std::vector<unsigned int> &pfs):_H(pfs),_hp(pfs.size(),0),_hv({pfs}),_x({1}){}
const Ham& operator++() {
for (int i=0; i<_H.size(); i++) for (;_hv[i]<=_x.back();_hv[i]=_x[++_hp[i]]*_H[i]);
_x.push_back(_hv[0]);
for (int i=1; i<_H.size(); i++) if (_hv[i]<_x.back()) _x.back()=_hv[i];
return *this;
}
};
int main() {
int count = 1;
for (unsigned int i : Ham({2,3,5})) {
if (count <= 62) std::cout << i << ' ';
if (count++ == 1691) {
std::cout << "\nThe one thousand six hundred and ninety first Hamming Number is " << i << std::endl;
break;
}
}
return 0;
}

复制代码

hujunhua · 发表于 2013-10-23 02:37:33

本帖最后由 hujunhua 于 2013-10-23 16:22 编辑

list={};(*预备一个空表来装要求的序列*)
X={2,3,5};(*倍乘因子表*)
checkb={1,1,1};(*备选框*)
pointer={0,0,0};(*倍乘因子指针表，指向已得的表list中的位置*)
Do[AppendTo[list,Min@checkb];(*将备选框的最小数追加到表list末*)
dp=1-Sign@(checkb-Last@list);(*用1标记上一步所选最小值在备选框中的位置，其它位置标为0，作为指针的移动表*)
pointer+=dp;(*按指针移动表移动指针*)
checkb+=(Thread[list[[pointer]]*dp*X]-Thread[checkb*dp]),(*将有移动的指针所指的表list中的项乘以对应的倍乘因子，以所得结果按标记更新备选框*)
{200}](*返回从备选框中选最小数入表的那一步，直至表长达到200*)
Print[list]

复制代码

是由楼上链结中的程序改编的，效率当然也很高喽。

mathematica · 发表于 2013-10-23 09:29:59

这个,每一项必然是前面几项中某一项乘以2或者3或者5而得到的

liangbch · 发表于 2013-10-23 17:10:06

这个题目的另一个说法是“求5-smooth数",百度知道(34#)给出的代码简单而有效，但有其局限性。用此代码求更大范围或者更大的质数的smooth数就不行了。其主要局限性是，
1.这个算法必须从1开始求所有smooth数，而不能求某一区间的smooth数。
2.这个算法所有求得的smooth数必须保存在内存。更大范围或者更大的质数的smooth数的个数会很多，可能导致内存不够用。不得不使用内存映射文件，或文件来存储,这会导致程序很复杂。

我正好需要计算大量的smooth数，鉴于磁盘和内存的限制，我不得不使用普通的做法，和素数的分段筛很类似，可分段求smooth数，不过求得需要的smooth数后需要排序操作，另外，我的这个算法好处是，便于分布式计算。

许多算法看上去很美，但很有局限性。和此算法一个类似的例子是《一个线性O(N)时间复杂度求素数筛法》(http://www.cnblogs.com/lzhenf/archive/2011/12/24/2300271.html),但使用后发现，虽然复杂度确实很小，但实际性能不理想，因为程序中对内存的访问局部性很差，导致cache命中率很低，限于这个算法特点，所有的所有素数必须保存在内存，故不能求很大范围内的素数。

下表是我2^64以内的p1到p32-smooth数的个数。

n	Pn	2^64以内的 Pn-smooth数的个数
1	2	64
2	3	1344
3	5	13282
4	7	85348
5	11	378554
6	13	1381206
7	17	4164316
8	19	11169655
9	23	26587554
10	29	56643654
11	31	113399428
12	37	210499995
13	41	370466869
14	43	627484292
15	47	1020181037
16	53	1591363196
17	59	2395982361
18	61	3527728592
19	67	5050259345
20	71	7077350886
21	73	9759368423
22	79	13183354386
23	83	17528499627
24	89	22918160232
25	97	29457339259
26	101	37424970530
27	103	47121505151
28	107	58735924658
29	109	72644250191
30	113	89058878497
31	127	107598623619
32	131	129070232683

mathematica · 发表于 2013-10-23 17:42:17

#include <stdio.h>
#include <limits.h>
#include <string.h>
#include <malloc.h>
#define MAX_LEVEL 32
typedef unsigned long DWORD;
#if defined(_MSC_VER)
typedef unsigned __int64 UINT64;
#define UINT64_FMT "%I64u"
#elif defined(__GNUC__)
typedef unsigned long long UINT64;
#define UINT64_FMT "%llu"
#else
#error don't support this compiler
#endif
const DWORD primes[]=
{
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
};
typedef struct _searchContent_tag
{
struct search_tag
{
UINT64 PI; // product of some numbers
UINT64 max;
int idx; // the index in array pArr
int maxIdx; // the index of last element
UINT64* pArr; // the address of primes power table
}arr[MAX_LEVEL];
UINT64 min;
UINT64 max;
int maxLevel; //level
}SRCH_CONTENT_ST;
void initSearchContent(int level, UINT64 min, UINT64 max, SRCH_CONTENT_ST* p)
{
UINT64 prime,limit,pp;
int i,c;
memset(p,0,sizeof(SRCH_CONTENT_ST));
p->maxLevel=level-1;
p->max=max; p->min=min;
p->arr[level-1].PI=1;
p->arr[level-1].max=max;
for (i=0;i<level;i++)
{
prime=primes[i];
limit=ULLONG_MAX/prime;
c=1;pp=1;
while (pp<=limit) //pp is power of prime
{ pp*=prime;c++; }
p->arr[i].pArr=(UINT64 *)malloc( sizeof(UINT64) * (c+1));
p->arr[i].pArr[0]=pp=1; c=1;
while (pp<=limit)
{ pp*=prime;p->arr[i].pArr[c++]=pp; }
p->arr[i].maxIdx= c-1;
}
}
void freeSearchContent( SRCH_CONTENT_ST* p)
{
int i;
for (i=0;i<=p->maxLevel;i++)
{
if ( p->arr[i].pArr!=NULL)
free(p->arr[i].pArr);
}
memset(p,0,sizeof(SRCH_CONTENT_ST));
}
//recursive form
void searchSmoothNums1( SRCH_CONTENT_ST* p, int order)
{
int idx;
/*
recurrence formula
p->arr[order].PI= 1;
p->arr[i-1].PI = arr[i].pArr[ arr[i].idx] * arr[i+1].PI;
p->arr[i].max= p->max / p->arr[i].PI;
*/
if (order==0)
{
UINT64 min;
UINT64 prod;
min= (p->min + p->arr[0].PI -1)/p->arr[0].PI; // min >= p->min/PI( when p->min % PI==0, min= p->min/PI)
if ( min > p->arr[0].max)
return;
idx=0;
while (idx < p->arr[0].maxIdx && p->arr[0].pArr[idx]<min )
idx++;
while ( p->arr[0].pArr[idx] <= p->arr[0].max && idx <= p->arr[0].maxIdx )
{
prod= p->arr[0].PI * p->arr[0].pArr[idx];
idx++;
printf(UINT64_FMT"\n",prod);
}
}
else
{
if ( p->arr[order].max==1 )
{
printf(UINT64_FMT"\n",p->arr[order].PI);
return ;
}
for (idx=0;
idx<= p->arr[order].maxIdx && p->arr[order].pArr[idx] <= p->arr[order].max;
idx++ )
{
p->arr[order].idx=idx;
p->arr[order-1].PI = p->arr[order].PI * p->arr[order].pArr[idx];
p->arr[order-1].max= p->max / p->arr[order-1].PI;
searchSmoothNums1(p,order-1); //call itselft
}
}
}
#define M_FILL 0 //fill mode
#define M_INC 1 //increase mode
// iterative form
// search all p[order+1] smooth number which between p->min and p->max
// For exmaple
// order=0, search 2-smooth numbers
// order=1, search 3-smooth numbers
// order=2, search 5-smooth number
void searchSmoothNums2( SRCH_CONTENT_ST* p, int order)
{
UINT64 prod;
UINT64 PI;
UINT64 min;
int idx;
int i=order;
int mode=M_FILL;
/*
recurrence formula
p->arr[order].PI= 1;
p->arr[i-1].PI = arr[i].pArr[ arr[i].idx] * arr[i+1].PI;
p->arr[i].max = p->max / arr[i].PI;
*/
p->arr[order].PI=1;
p->arr[order].max= p->max;
while (i<=order)
{
if (mode==M_FILL) //fill mode
{
PI = p->arr[i].PI;
if (i>0)
{
p->arr[i].idx=0;
p->arr[i-1].PI = p->arr[i].PI; //p->arr[i-1].PI= p->arr[i].PI;
p->arr[i-1].max= p->arr[i].max;
i--; //forward to next level
}
else // i==0
{
min= (p->min + PI-1)/PI; // min >= p->min/PI( when p->min % PI==0, min= p->min/PI)
if ( min > p->arr[0].max)
{
i++; // backtrace
mode=M_INC;
}
else
{
idx=0;
while (idx < p->arr[0].maxIdx && p->arr[0].pArr[idx]<min )
idx++;
if ( p->arr[0].pArr[idx] > p->arr[0].max || idx > p->arr[0].maxIdx)
i++; // backtrace
else //in this case, varible i don't change
{
p->arr[0].idx=idx;
printf(UINT64_FMT"\n",PI * p->arr[0].pArr[idx]);
}
mode=M_INC; //change to increase mode
}
}
}
else //is increase mode
{
idx= p->arr[i].idx + 1;
if (p->arr[i].pArr[idx] > p->arr[i].max || idx> p->arr[i].maxIdx)
i++; // backtrace, mode don't change, still is increase mode
else
{
p->arr[i].idx=idx;
prod= p->arr[i].PI * p->arr[i].pArr[idx];
if (i==0)
printf(UINT64_FMT"\n",prod);
else
{
p->arr[i-1].max= p->max / prod;
if ( p->arr[i-1].max ==1)
{
printf(UINT64_FMT"\n",prod);
i++; // backtrace, mode don't change, still is increase mode
}
else
{
p->arr[i-1].PI=prod;
mode=M_FILL; // change mode to fill mode
i--; // forward to next level
}
}
}
}
}
}
void ref_searchSmoothNums(int order,UINT64 min, UINT64 max)
{
UINT64 i,m;
int j,isSmooth;
for (i=min;i<=max;i++)
{
for (m=i,isSmooth=0,j=0;j<order;j++)
{
while (m % primes[j]==0) m/=primes[j];
if (m==1)
{ isSmooth=1; break; }
}
if (isSmooth)
printf(UINT64_FMT"\n",i);
}
}
void test()
{
SRCH_CONTENT_ST searchContent;
UINT64 nMin,nMax;
int cs;//case
for (cs=0;cs<=1;cs++)
{
if (cs==0)
{ nMin=1; nMax=100; }
else if (cs==1)
{ nMin=100; nMax=200; }
initSearchContent(4,nMin,nMax,&searchContent);
printf("\n\nSearch smooth number ref_searchSmoothNums-------\n");
ref_searchSmoothNums(4,nMin,nMax);
printf("\n\nSearch smooth number searchSmoothNums1-------\n");
searchSmoothNums1( &searchContent,4-1);
printf("\n\nSearch smooth number searchSmoothNums2-------\n");
searchSmoothNums2( &searchContent,4-1);
freeSearchContent( &searchContent);
}
}
int main(int argc, char* argv[])
{
test();
return 0;
}

复制代码

1.什么是Smooth 数

如果一个整数的所有素因子都不大于B，我们称这个整数为B-Smooth数。用符号Pn表示表示第n个素数，例P1=2,P2=3,P3=5. 我们称素因子不大于Pn的Smooth数叫Pn-Smooth数。关于Smooth数更多的信息，请请参阅 Wiki词条 Smooth number http://en.wikipedia.org/wiki/Smooth_number.  P1 到 P9-Smooth 数可参阅下面的数列，这些数列在On-Line Encycolpedia of Integer Sequences http://oeis.org 给出

 2-smooth numbers: A000079 (2^i, i>=0)

 3-smooth numbers: A003586 (2^i ×3^j,i,j>=0)

 5-smooth numbers: A051037 (2^i ×3^j×5^k, i,j,k>=0)

 7-smooth numbers: A002473 (2^i×3^j×5^k×7^l, i,j,k,l>=0)

 11-smooth numbers: A051038 (etc...)

 13-smooth numbers: A080197

 17-smooth numbers: A080681

 19-smooth numbers: A080682

 23-smooth numbers: A080683

2. 如何得到smooth数

  基本上，我们有2种方法来生成一定范围内的smooth数。

1. 穷举法。基本方法是对指定范围内的每个数做分解素因数，如果其所有素因子都小于等于Pn,那个这个数就是Pn-smooth数。这种方法的优点是简单，缺点是速度较慢，只适合检查小范围内的整数。

2．构造法。其基本方法是用一个数组表示各个素因子的指数。例arr是一个数组，arr[1]表示素数p1（2）的指数，arr[2]表示素数p2（3）的指数，arr[3]表示素数p3（5）的指数. 欲求min和max之间的所有pn-smooth数, 则需要列举枚举这个数组的各个元素的取值。数组的各个元素的值需满足如下2个条件。

1. p[ i]>=0, 且p[ i] ^ arr[ i] <= max.  （p[ i]表示第i个素数，i>=1,i<=n）

2. min <=p[1]^arr[1] * p[2]^arr[2] * p[3]^arr[3] * …p[n]*arr[n] <=  max

对于一个包含n个元素的数组arr, 每个元素arr[ i] （1<=i<=n）各取一个数，就得到一个n维向量。如arr[1]=0, arr[2]=0, arr[3]=0可用3维向量(0,0,0)来表示，而arr[1]=1, arr[2]=0, arr[3]=0 则用3维向量(1,0,0)来表示。为了使得取得的所有n维向量既不重复也不遗漏。我们这里定义了比较2个n维向量的大小的规则

1.  如果2个n维向量a,b的所有分量都是相同的，我们说这两个n维向量是相同的。

2.  对于2个n维向量a,b（a的各个分量为a[1],a[2], …, a[n]）,如果a[n]>b[n]我们说a>b.

3.  如果存在某个i，使得（1）对于所有的j (i<j<=n)都有a[j]=b[j]. (2) a[ i]>b[ i],那么我们说n维向量a>b.

  我们列举每个n维向量时，总是从小到到的进行，先是(0,0,0), 然后是(1,0,0), 然后是(2,0,0)…，然后是(0,1,0), (1,1,0), (2,1,0).

未完待续

http://blog.csdn.net/liangbch/article/details/12186099

chyanog · 发表于 2013-10-23 17:48:53

本帖最后由 chyanog 于 2013-10-23 17:54 编辑

hujunhua 发表于 2013-10-23 02:37
是由楼上链结中的程序改编的，效率当然也很高喽。

当list很大的时候，AppendTo会越来越慢，不如用Internal`Bag

Module[{list, X, dp, checkb, pointer},
list = {};
X = {2, 3, 5};
checkb = {1, 1, 1};
pointer = {0, 0, 0};
Do[AppendTo[list, Min@checkb];
dp = 1 - Sign@(checkb - Last@list);
pointer += dp;
checkb += (list[[pointer]]*dp*X - checkb*dp)
, {2 10^4}];
list
] // Hash // AbsoluteTiming

复制代码

{5.515315, 1434246505}

\$ContextPath = Union@Append[\$ContextPath, "Internal`"];
Module[{bag, X, dp, checkb, pointer},
bag = Bag[];
X = {2, 3, 5};
checkb = {1, 1, 1};
pointer = {0, 0, 0};
Do[
StuffBag[bag, Min@checkb];
dp = 1 - Sign@(checkb - BagPart[bag, -1]);
pointer += dp;
checkb += (BagPart[bag, pointer]*dp*X - checkb*dp),
{2 10^4}];
bag~BagPart~All
] // Hash // AbsoluteTiming

复制代码

{0.905052, 1434246505}

账号		自动登录	找回密码
密码			欢迎注册

[擂台] 求形如 2^r3^s5^t 的整数序列

点评

点评

点评

Mathematica程序，不显含 if 语句。

点评

点评

点评

[擂台] 求形如 2^r*3^s*5^t 的整数序列

点评

点评

点评

Mathematica程序，不显含 if 语句。

点评

点评

点评

[擂台] 求形如 2^r3^s5^t 的整数序列