Mathematica 求幂并查找指定系数

发布于 2024-11-17 03:46:04 字数 459 浏览 5 评论 0原文

我有以下代码，它完全按照我想要的方式执行，只是速度慢得离谱。我不会那么烦恼，除了当我“手动”处理代码时，即我将其分成几部分并单独执行它们时，它几乎是瞬时的。

这是我的代码：

Coefficient[Product[Sum[x^(j*Prime[i]), {j, 0, Floor[q/Prime[i]]}], 
                        {i, 1, PrimePi[q]}], x, q]

为了清楚起见添加图片：

在此处输入图像描述

我认为它正在尝试优化总和，但我没有把握。有办法阻止吗？

另外，由于我所有的系数都是正数，并且我只想要第 x^q 个系数，有没有办法让 Mathematica 丢弃所有大于该指数的指数，而不与这些指数进行所有乘法？

原文

I have the following code, and it does exactly what I want it to do, except that it is ridiculously slow. I would not be so bothered, except that when I process the code "manually", i.e., I break it into parts and do them individually, it's near instantaneous.

Here is my code:

Coefficient[Product[Sum[x^(j*Prime[i]), {j, 0, Floor[q/Prime[i]]}], 
                        {i, 1, PrimePi[q]}], x, q]

Picture added for clarity:

enter image description here

I think it is trying to optimize the sum, but am not sure. Is there a way to stop that?

In addition, since all my coefficients are positive, and I only want the x^qth one, is there a way to get Mathematica to discard all exponents that are larger than that and not do all the multiplication with those?

分享到QQ

分享到微博

如果你对这篇内容有疑问，欢迎到本站社区发帖提问参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

发布评论

需要登录才能够评论，你可以免费注册一个本站的账号。

星星的轨迹 2024-11-24 03:46:04

我可能会误解您想要的内容，但是由于系数取决于 q，我假设您希望针对特定的 q 对其进行评估。由于我（像您一样）怀疑优化乘积和总和需要花费时间，因此我重写了它。你有类似的东西：

With[{q = 80}, Coefficient[\!\(
\*UnderoverscriptBox[\(\[Product]\), \(i = 1\), \(PrimePi[q]\)]\((
\*UnderoverscriptBox[\(\[Sum]\), \(j = 0\), \(\[LeftFloor]
\*FractionBox[\(q\), \(Prime[i]\)]\[RightFloor]\)]
\*SuperscriptBox[\(x\), \(j*Prime[i]\)])\)\), x, q]] // Timing
(*
-> {8.36181, 10003}
*)

我用纯粹的结构操作重写为

With[{q = 80},
 Coefficient[Times @@ 
 Table[Plus @@ Table[x^(j*Prime[i]), {j, 0, Floor[q/Prime[i]]}],
        {i, 1, PrimePi[q]}], x, q]] // Timing
(*
-> {8.36357, 10003}
*)

（这只是建立了一个术语列表，然后将它们相乘，因此不执行符号分析）。

建立多项式是瞬时的，但它有几千个项，因此可能发生的情况是 Coefficient 花费大量时间来确保它具有正确的系数。实际上，您可以通过展开多项式来解决这个问题。因此：

 With[{q = 80}, Coefficient[Expand[\!\(
 \*UnderoverscriptBox[\(\[Product]\), \(i = 1\), \(PrimePi[q]\)]\((
 \*UnderoverscriptBox[\(\[Sum]\), \(j = 0\), \(\[LeftFloor]
 \*FractionBox[\(q\), \(Prime[i]\)]\[RightFloor]\)]
 \*SuperscriptBox[\(x\), \(j*Prime[i]\)])\)\)], x, q]] // Timing
 (*
 -> {0.240862, 10003}
 *)

它也适用于我的方法。

总而言之，只需将 Expand 粘贴在表达式前面和获取系数之前即可。

I may be misunderstanding what you want but, as the coefficient will depend on q, I assume you want it evaluated for specific q. Since I suspected (like you) that the time is taken to optimise the produt and sum, I rewrote it. You had something like:

With[{q = 80}, Coefficient[\!\(
\*UnderoverscriptBox[\(\[Product]\), \(i = 1\), \(PrimePi[q]\)]\((
\*UnderoverscriptBox[\(\[Sum]\), \(j = 0\), \(\[LeftFloor]
\*FractionBox[\(q\), \(Prime[i]\)]\[RightFloor]\)]
\*SuperscriptBox[\(x\), \(j*Prime[i]\)])\)\), x, q]] // Timing
(*
-> {8.36181, 10003}
*)

which I rewrote with purely structural operations as

With[{q = 80},
 Coefficient[Times @@ 
 Table[Plus @@ Table[x^(j*Prime[i]), {j, 0, Floor[q/Prime[i]]}],
        {i, 1, PrimePi[q]}], x, q]] // Timing
(*
-> {8.36357, 10003}
*)

(this just builds up a list of the terms and then multiplies them, so no symbolic analysis is performed).

Just building up the polynomial is instantaneous, but it has a few thousand terms, so what is probably happening is that Coefficient spends a lot of time to make sure it has the right coefficient. Actually you can solve this by Expanding the polynomial. Thus:

 With[{q = 80}, Coefficient[Expand[\!\(
 \*UnderoverscriptBox[\(\[Product]\), \(i = 1\), \(PrimePi[q]\)]\((
 \*UnderoverscriptBox[\(\[Sum]\), \(j = 0\), \(\[LeftFloor]
 \*FractionBox[\(q\), \(Prime[i]\)]\[RightFloor]\)]
 \*SuperscriptBox[\(x\), \(j*Prime[i]\)])\)\)], x, q]] // Timing
 (*
 -> {0.240862, 10003}
 *)

and it also works for my method.

So to summarise, just stick Expand in front of the expression and before you take the coefficient.

回复收藏 0 原文

兮子 2024-11-24 03:46:04

我认为原始代码速度慢的原因是因为 Coefficient 甚至可以使用非常大的表达式 - 如果天真地扩展，这些表达式将无法装入内存。

这是原始多项式：

poly[q_, x_] := Product[Sum[ x^(j*Prime[i]), 
                            {j, 0, Floor[q/Prime[i]]}], {i, 1, PrimePi[q]}]

看看对于不太大的 q，扩展多项式会占用更多内存并且变得相当慢：

In[2]:= Through[{LeafCount, ByteCount}[poly[300, x]]] // Timing
        Through[{LeafCount, ByteCount}[Expand@poly[300, x]]] // Timing
Out[2]= { 0.01, { 1859,   55864}}
Out[3]= {25.27, {77368, 3175840}}

现在让我们以 3 种不同的方式定义系数并为它们计时

coeff[q_]    := Module[{x}, Coefficient[poly[q, x], x, q]]
exCoeff[q_]  := Module[{x}, Coefficient[Expand@poly[q, x], x, q]]
serCoeff[q_] := Module[{x}, SeriesCoefficient[poly[q, x], {x, 0, q}]]

In[7]:= Table[   coeff[q],{q,1,30}]//Timing
        Table[ exCoeff[q],{q,1,30}]//Timing
        Table[serCoeff[q],{q,1,30}]//Timing
Out[7]= {0.37,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}
Out[8]= {0.12,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}
Out[9]= {0.06,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}

In[10]:=    coeff[100]//Timing
          exCoeff[100]//Timing
         serCoeff[100]//Timing
Out[10]= {56.28,40899}
Out[11]= { 0.84,40899}
Out[12]= { 0.06,40899}

所以 SeriesCoefficient 绝对是正确的选择。当然除非你是
比我更擅长组合数学，并且您知道以下素数分配公式
(oeis)

In[13]:= CoefficientList[Series[1/Product[1-x^Prime[i],{i,1,30}],{x,0,30}],x]
Out[13]= {1,0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}

In[14]:= f[n_]:=Length@IntegerPartitions[n,All,Prime@Range@PrimePi@n]; Array[f,30]
Out[14]= {0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}

I think that the reason that the original code is slow is because Coefficient is made to work even with very large expressions - ones that would not fit into the memory if naively expanded.

Here's the original polynomial:

poly[q_, x_] := Product[Sum[ x^(j*Prime[i]), 
                            {j, 0, Floor[q/Prime[i]]}], {i, 1, PrimePi[q]}]

See how for not too large q, expanding the polynomial takes up a lot more memory and becomes fairly slow:

In[2]:= Through[{LeafCount, ByteCount}[poly[300, x]]] // Timing
        Through[{LeafCount, ByteCount}[Expand@poly[300, x]]] // Timing
Out[2]= { 0.01, { 1859,   55864}}
Out[3]= {25.27, {77368, 3175840}}

Now let's define the coefficient in 3 different ways and time them

coeff[q_]    := Module[{x}, Coefficient[poly[q, x], x, q]]
exCoeff[q_]  := Module[{x}, Coefficient[Expand@poly[q, x], x, q]]
serCoeff[q_] := Module[{x}, SeriesCoefficient[poly[q, x], {x, 0, q}]]

In[7]:= Table[   coeff[q],{q,1,30}]//Timing
        Table[ exCoeff[q],{q,1,30}]//Timing
        Table[serCoeff[q],{q,1,30}]//Timing
Out[7]= {0.37,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}
Out[8]= {0.12,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}
Out[9]= {0.06,{0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}}

In[10]:=    coeff[100]//Timing
          exCoeff[100]//Timing
         serCoeff[100]//Timing
Out[10]= {56.28,40899}
Out[11]= { 0.84,40899}
Out[12]= { 0.06,40899}

So SeriesCoefficient is definitely the way to go. Unless of course you're
a bit better at combinatorics than me and you know the following prime partition formulae
(oeis)

In[13]:= CoefficientList[Series[1/Product[1-x^Prime[i],{i,1,30}],{x,0,30}],x]
Out[13]= {1,0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}

In[14]:= f[n_]:=Length@IntegerPartitions[n,All,Prime@Range@PrimePi@n]; Array[f,30]
Out[14]= {0,1,1,1,2,2,3,3,4,5,6,7,9,10,12,14,17,19,23,26,30,35,40,46,52,60,67,77,87,98}

回复收藏 0 原文

~没有更多了~