【问题标题】:How to make code for sum of 2 integrals in Mathematica如何在 Mathematica 中为 2 个积分之和编写代码
【发布时间】:2012-01-11 09:31:05
【问题描述】:

我正在尝试快速解决以下问题:

f[r_] := Sum[(((-1)^n (2 r - 2 n - 7)!!)/(2^n n! (r - 2 n - 1)!))
             * x^(r - 2*n - 1), 
         {n, 0, r/2}]; 

Nw := Transpose[Table[f[j], {i, 1}, {j, 5, 200, 1}]]; 

X1 = Integrate[Nw . Transpose[Nw], {x, -1, 1}]

我可以通过这段代码快速得到答案:

$starttime = AbsoluteTime[]; Quiet[LaunchKernels[]]; 
DIM = 50; 
Print["$Version = ", $Version, "  |||  ", 
      "Number of Kernels : ", Length[Kernels[]]]; 

Nw = Transpose[Table[f[j], {i, 1}, {j, 5, DIM, 1}]]; 
nw2 = Nw.Transpose[Nw]; 
Round[First[AbsoluteTiming[nw3 = ParallelMap[Expand, nw2]; ]]] 

intrule = (pol_Plus)?(PolynomialQ[#1, x]&) :> 
      (Select[pol, !FreeQ[#1, x] & ] /. 
         x^(n_.) /; n > -1 :> ((-1)^n + 1)/(n + 1)) + 2*(pol /. x -> 0)]); 

Round[First[AbsoluteTiming[X1 = ParallelTable[row /. intrule, {row, nw3}]; ]]]

X1 
Print["overall time needed in seconds: ", Round[AbsoluteTime[] - $starttime]]; 

但是如果我需要解决以下问题,其中 a 和 b 是已知常量,我该如何管理这段代码?

       X1 = a Integrate[Nw.Transpose[Nw], {x, -1, 0.235}]
          + b Integrate[Nw.Transpose[Nw], {x,  0.235,1}]; 

【问题讨论】:

  • 这个问题不是特别清楚,格式也很差,让人难以阅读,因此不鼓励人们帮助你。 you 可以尝试修复它吗(这包括代码周围的标题和文本)? -- 请注意,好的问题会鼓励好的答案,然后每个人都会受益!
  • 相关问题:SO/8021501 & SU/315337
  • 西蒙有什么不清楚的地方?我首先需要 X1 = Integrate[Nw 。转置[Nw], {x, -1, 1}];现在我需要 X1 = aIntegrate[Nw 。 Transpose[Nw], {x, -1, 0.235}]+bIntegrate[Nw .Transpose[Nw], {x, 0.235,1}];所以尝试用给出的代码来获得这两个积分,这就是全部
  • 可以理解,只是不清楚。
  • 欢迎回来。我有一些问题。您介意我进行一些格式化,包括删除Print 语句吗?我想只包含核心代码,因为它有助于查看代码在做什么。其次,在Nw 的定义中,您使用:=,它会在每次调用Nw 时重新评估。这是故意的,为什么?第三,Nw中的Table只生成一行(i总是1),那么你有没有i的值不止一个的情况?

标签: wolfram-mathematica integration rotational-matrices


【解决方案1】:

这是一个对多项式进行定积分的简单函数

polyIntegrate[expr_List, {x_, x0_, x1_}] := polyIntegrate[#, {x, x0, x1}]&/@expr
polyIntegrate[expr_, {x_, x0_, x1_}] := Check[Total[# 
  Table[(x1^(1 + n) - x0^(1 + n))/(1 + n), {n, 0, Length[#] - 1}]
  ]&[CoefficientList[expr, x]], $Failed, {General::poly}]

就其适用范围而言,这比使用Integrate 快大约100 倍。对于您的问题,这应该足够快。如果不是,则可以并行化。

f[r_] := Sum[(((-1)^n*(2*r - 2*n - 7)!!)/(2^n*n!*(r - 2*n - 1)!))*
   x^(r - 2*n - 1), {n, 0, r/2}];
Nw = Transpose[Table[f[j], {i, 1}, {j, 5, 50, 1}]];

a*polyIntegrate[Nw.Transpose[Nw], {x, -1, 0.235}] + 
   b*polyIntegrate[Nw.Transpose[Nw], {x, 0.235, 1}] // Timing // Short

(* Returns: {7.9405,{{0.0097638 a+0.00293462 b,<<44>>,
             -0.0000123978 a+0.0000123978 b},<<44>>,{<<1>>}}} *)

【讨论】:

  • 我正要告诫你没有为polyIntegrate 提供足够的保护来处理expr == constantexpr == 1/x,但显然polyIntegrate 处理第一个很好,第二个处理得很好从来没有出现过,尽管看起来会。顺便说一句,+1。
  • 谢谢你亲爱的西蒙,无论如何如何并行化?
  • @Simon 我刚刚使用了旧的并行化代码,但在前面使用 poly 对其进行了改进,并且运行速度非常快。非常感谢!
  • @rcollyer: 如果expr 不是x 中的多项式,CoefficientList 会报错。我现在为此添加了Check,因此polyIntegrate 在给定非多项式时将返回$Failed。 (昨天懒得加了!)
  • @Simon,我注意到... :P 或者,您可以使用PolynomialQ[#,x]&amp;
猜你喜欢
  • 1970-01-01
  • 2019-04-30
  • 1970-01-01
  • 1970-01-01
  • 2016-02-19
  • 1970-01-01
  • 2019-07-03
  • 1970-01-01
  • 2017-04-25
相关资源
最近更新 更多