2k^2+2k+1是完全平方数？

waynebuaa

来一个简单的题玩玩。


如果k是正整数，且2k^2+2k+1是完全平方数，这样的k有多少？

看谁的算法好，看谁的代码少，看谁找的多

jimogsh

原来不知道一个数是不是完全平方数不知道有没有现成的函数，后来发现了一个MathchQ，就现学现用了一下，不算简洁，也肯定没有找完，先拿出来晒一下，呵呵

1. In[12]:= f[k_] := 2 k^2 + 2 k + 1; For[k = 1, k <= 1000000000, k++,
   
2. If[MatchQ[Sqrt[f[k]], _Integer], Print[k]]]
   
3. 
   
4. During evaluation of In[12]:= 3
   
5. 
   
6. During evaluation of In[12]:= 20
   
7. 
   
8. During evaluation of In[12]:= 119
   
9. 
   
10. During evaluation of In[12]:= 696
    
11. 
    
12. During evaluation of In[12]:= 4059
    
13. 
    
14. During evaluation of In[12]:= 23660
    
15. 
    
16. During evaluation of In[12]:= 137903
    
17. 
    
18. During evaluation of In[12]:= 803760
    
19. 
    
20. Out[12]= $Aborted

复制代码

waynebuaa

一个一个的检验一个数是否能开平方，速度估计是最慢的了
Floor[(Sqrt[2] + 1)^(2 n - 1)/4] /. n -> Range[50]

waynebuaa

我最开始也是用你这种方法，后来换成了Reduce，再后来用上了FindLinearRecurrence，感觉把这个问题搞的比较明白清楚了，可在数学论坛里请教了人家之后，才发现还有几种迅猛的算法，:victory:

FreddyMusic

从数学上讲，还是个无穷数列，没看懂中间如何推导的，你把推导过程说清楚。

不明白为何不用 Map 而喜欢 ReplaceAll，我的代码比你的少几个字符。

1. 
   
2. 
   
3. Floor[(Sqrt[2] + 1)^(2 # - 1)/4] & /@ Range[50]
   
4. 
   

复制代码

waynebuaa

这样闹没意义，......

硬要这么闹也轮不到你给所谓的Map

1. Floor[(Sqrt[2] + 1)^(2 Range[50] - 1)/4]

复制代码

jimogsh

也许可以通过找交集的方法来“找”，当然这比找到一个通解麻烦多了，不过我想还是比较快的

1. In[2]:= f[k_] := 2 k^2 + 2 k + 1; NSolve[f[x] == 1000000, x]
   
2. 
   
3. Out[2]= {{x -> -707.607}, {x -> 706.607}}
   
4. 
   
5. In[6]:= lis1 = Table[k^2, {k, 0, 706}]; lis2 =
   
6. Table[f[k], {k, 0, 1000}]; Intersection[lis1, lis2]
   
7. 
   
8. Out[6]= {1, 25, 841, 28561}

复制代码

当然这里找出的是完全平方数k^2而非k本身

FreddyMusic

way, you are not exactly see the future and timing.
Map is fundmental and it's good, later you can move to ParallelMap with parallel kernels.
For example: when it goes to big number, the advantage shows.
Can your code do this ?

1. 
   
2. 
   
3. In[289]:= ClearSystemCache[];
   
4. Floor[(Sqrt[2]+1)^(2Range[5000]-1)/4];//AbsoluteTiming
   
5. Out[290]= {5.3906250,Null}
   
6. In[291]:= ClearSystemCache[];
   
7. Map[Floor[(Sqrt[2]+1)^(2 #-1)/4]&,Range[5000]];//AbsoluteTiming
   
8. Out[292]= {5.4687500,Null}
   
9. In[293]:= ClearSystemCache[];
   
10. ParallelMap[Floor[(Sqrt[2]+1)^(2 #-1)/4]&,Range[5000]];//AbsoluteTiming
    
11. Out[294]= {2.0781250,Null}
    
12. 
    

复制代码

waynebuaa

我的机子是七年前的，玩不起Parallel computing
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

simwe论坛速度好慢啊，我都没耐性了。