加法数论

kastin

对于不定方程$x_1^k+x_2^k+x_3^k+...+x_s^k=m$，以及正整数k,s,m，求：
1.方程的正整数解的个数；
2.方程的非负整数解的个数；
3.方程偏序解的个数(即满足$x_1<x_2<...<x_s$的解)。

如果能给出公式请给出（不一定是封闭形式，可以使递归式），若不能，请给出高效算法，并给出例子验证。

kastin

假设问题(1)的解的个数为$A^+(,k,s,m)$, (2)的解的个数为$A^N(k,s,m)$, (3)的解的个数为$A^o(k,s,m)$，那么有
a) $A^N(k,s,m)=A^+(k,1,m)+A^+(k,2,m)+...+A^+(k,s,m)$
因此只要求出$A^+(k,s,m)$即可.
对于$A^+(k,s,m)$，有递归关系
b) $A^+(k,s,m)=A^+(k,s-1,m-1)+A^+(k,s-1,m-2^k)+cdots+A^+(k,s-1,m-\[m^{1/k}\]^k)+\[m^{1/k}\]^k$,  其中$\[·\]$表示高斯取整函数。
考虑到$A^+(k,s,m)=0$, 如果$m<s$, 因此b)中项数并不全部非零。因此对右端每一项反复使用上式，使得s递降，那么可以得到$A^+(k,s,m)$了。不过这种方法对于于s,m不大时候，使用计算机编程可以求解。当s,k或m过大，递归次数剧增使得计算机开销过载。因此这个方法不是最好的办法。

kastin

问题2的公式有了，但是无法用mathematica实现。因为中间有个算子运算。

这里的$\Theta_k(t)=(-1)^t(C_t^{1^k}+C_t^{2^k}+cdots+C_t^{r^k}).$  这里$r=[t^(1/k)].

kastin

第三问有人有思路吗？
我有种思路，只是更加复杂了：
由于$x_1<x_2<cdots<x_s$，因此进行下面的变换，可得等价形式
$X_1=x_1, X_2=x_1+x_2, X3=x_1+x2+x_3, cdots, X_s=x_1+x_2+,cdots,+x_s$
将$X_1,..,X_s$代入原方程，那么$(x_1,x_2,cdots,x_s)$正整数解的个数即为问题(3)中解的个数。

kastin

@wayne
我知道有啊，但是它用的是搜索算法，这个没意义啊。

wayne

kastin 发表于 2013-10-30 10:19
@wayne
我知道有啊，但是它用的是搜索算法，这个没意义啊。

其实动态规划已经足够快了

kastin

wayne 发表于 2013-10-30 17:36
其实动态规划已经足够快了

能否给个程序呢？mathematica的算法我不知道是什么，但是对于有些特殊输入，PowersRepresentations运行比较慢。

wayne

kastin 发表于 2013-11-1 18:56
能否给个程序呢？mathematica的算法我不知道是什么，但是对于有些特殊输入，PowersRepresentations运行比 ...

这个问题 我之前在论坛问过。
关于整数的等幂和拆分的计数问题
从理论上分析 ，基本上就你说的这些,我也没啥好点子 。
在程序实现上，我们可以将小的结果存储起来，以空间换时间(动态规划)，避免深度迭代和不必要的重复计算。
也可以这么考虑，我们不妨把所有的小于m的k次幂 存起来，解此线性的方程$\sum_{i=0}^{n}a_i*x_i =m$ ，其中$n^k<=m$，系数$a_i=i^k$,$x_i$之和为s

wayne

但是对于有些特殊输入，PowersRepresentations运行比较慢。

这种输入基本上都是 s比较大，同时k很小(2或者3)，因为此时 解的个数 本身就很多。

而 当k很大或者s很小的时候，该函数还是很快就能返回结果的

wayne

我试着写了个程序，计算10000表示成15个正整数的立方之和的解的个数，跟Mathematica自带的函数PowersRepresentations 做了对比，
我的代码花了58秒钟，返回8292个解，而PowersRepresentations只需40秒钟就返回了。
然后我搜索了下Mathematica的系统文件，发现PowersRepresentations函数被编译成mx这种二进制文件了，也就说，算法基本上是一致的。

1. Select[PowersRepresentations[10000, 15, 3], FreeQ[#, 0] &]

复制代码

1. m = 10000; s = 15; k = 3; r = Floor[m^(1/k)];
   
2. x /@ Range[r] /.
   
3.    Solve[Total[(Range[r]^3).(x /@ Range[r])] == m &&
   
4.      Total[x /@ Range[r]] == s && ((0 <= # <= r & /@ (x /@ Range[r])) /. List -> And),
   
5.     x /@ Range[r], Integers]
   

复制代码