2

Mathematica 7または8で積分を行う最良の方法は何ですか

NIntegrate[Exp[-x]/Sin[Pi x], {x, 0, 50}]

すべての整数に極があり、コーシー原理の値が必要です。アイデアは、0 から無限大までの積分の適切な近似を得ることです。

IntegrateのオプションがありますPrincipleValue -> True

I を使用NIntegrateすると、オプションExclusions -> (Sin[Pi x] == 0)を指定するか、手動で極を指定できます

NIntegrate[Exp[-x]/Sin[Pi x], Evaluate[{x, 0, Sequence@@Range[50], 50}]]

元のコマンドと上記の 2 つのNIntegrateトリックにより、結果が得られ60980 +/- 10ます。しかし、それらはすべてエラーを吐き出します。Mathematica がエラーを出さないように、この積分の信頼できる結果を迅速に得る最善の方法は何ですか?

4

3 に答える 3

7

サイモン、あなたの積分が収束していると信じる理由はありますか?

In[52]:= f[k_Integer, eps_Real] := 
 NIntegrate[Exp[-x]/Sin[Pi x], {x, k + eps, k + 1 - eps}]

In[53]:= Sum[f[k, 1.0*10^-4], {k, 0, 50}]

Out[53]= 2.72613

In[54]:= Sum[f[k, 1.0*10^-5], {k, 0, 50}]

Out[54]= 3.45906

In[55]:= Sum[f[k, 1.0*10^-6], {k, 0, 50}]

Out[55]= 4.19199

問題はx==0にあるようです。kの整数値の被積分関数k+epsをk+1-epsに分割します。

In[65]:= int = 
 Sum[(-1)^k Exp[-k ], {k, 0, Infinity}] Integrate[
   Exp[-x]/Sin[Pi x], {x, eps, 1 - eps}, Assumptions -> 0 < eps < 1/2]

Out[65]= (1/((1 + 
   E) (I + \[Pi])))E (2 E^(-1 + eps - I eps \[Pi])
     Hypergeometric2F1[1, (I + \[Pi])/(2 \[Pi]), 3/2 + I/(2 \[Pi]), 
     E^(-2 I eps \[Pi])] + 
   2 E^(I eps (I + \[Pi]))
     Hypergeometric2F1[1, (I + \[Pi])/(2 \[Pi]), 3/2 + I/(2 \[Pi]), 
     E^(2 I eps \[Pi])])

In[73]:= N[int /. eps -> 10^-6, 20]

Out[73]= 4.1919897038160855098 + 0.*10^-20 I

In[74]:= N[int /. eps -> 10^-4, 20]

Out[74]= 2.7261330651934049862 + 0.*10^-20 I

In[75]:= N[int /. eps -> 10^-5, 20]

Out[75]= 3.4590554287709991277 + 0.*10^-20 I

ご覧のとおり、対数特異点があります。

In[79]:= ser = 
 Assuming[0 < eps < 1/32, FullSimplify[Series[int, {eps, 0, 1}]]]

Out[79]= SeriesData[eps, 0, {(I*(-1 + E)*Pi - 
     2*(1 + E)*HarmonicNumber[-(-I + Pi)/(2*Pi)] + 
          Log[1/(4*eps^2*Pi^2)] - 2*E*Log[2*eps*Pi])/(2*(1 + E)*Pi), 
     (-1 + E)/((1 + E)*Pi)}, 0, 2, 1]

In[80]:= Normal[
  ser] /. {{eps -> 1.*^-6}, {eps -> 0.00001}, {eps -> 0.0001}}

Out[80]= {4.191989703816426 - 7.603403526913691*^-17*I, 
 3.459055428805136 - 
     7.603403526913691*^-17*I, 
 2.726133068607085 - 7.603403526913691*^-17*I}

上記のコードのEDITOut [79]は、eps-> 0の級数展開を示し、これら2つの対数項を組み合わせると、次のようになります。

In[7]:= ser = SeriesData[eps, 0, 
       {(I*(-1 + E)*Pi - 2*(1 + E)*HarmonicNumber[-(-I + Pi)/(2*Pi)] + 
              Log[1/(4*eps^2*Pi^2)] - 2*E*Log[2*eps*Pi])/(2*(1 + E)*
       Pi), 
         (-1 + E)/((1 + E)*Pi)}, 0, 2, 1]; 

In[8]:= Collect[Normal[PowerExpand //@ (ser + O[eps])], 
 Log[eps], FullSimplify]

Out[8]= -(Log[eps]/\[Pi]) + (
 I (-1 + E) \[Pi] - 
  2 (1 + E) (HarmonicNumber[-((-I + \[Pi])/(2 \[Pi]))] + 
     Log[2 \[Pi]]))/(2 (1 + E) \[Pi])

明らかに、-Log [eps]/Piはx==0の極から来ました。したがって、これを差し引くと、主値法が他の極に対してこれを行うのと同じように、最終的に有限値になります。

In[9]:= % /. Log[eps] -> 0

Out[9]= (I (-1 + E) \[Pi] - 
 2 (1 + E) (HarmonicNumber[-((-I + \[Pi])/(2 \[Pi]))] + 
    Log[2 \[Pi]]))/(2 (1 + E) \[Pi])

In[10]:= N[%, 20]

Out[10]= -0.20562403655659928968 + 0.*10^-21 I

もちろん、この結果を数値で検証することは困難ですが、私があなたの問題について行っていることをもっと知っているかもしれません。

編集2

この編集は、元の正則化された積分を計算するIn[65]入力を正当化するためのものです。私たちはコンピューティングしています

Sum[ Integrate[ Exp[-x]/Sin[Pi*x], {x, k+eps, k+1-eps}], {k, 0, Infinity}] ==  
  Sum[ Integrate[ Exp[-x-k]/Sin[Pi*(k+x)], {x, eps, 1-eps}], {k, 0, Infinity}] ==
  Sum[ (-1)^k*Exp[-k]*Integrate[ Exp[-x]/Sin[Pi*x], {x, eps, 1-eps}], 
       {k, 0, Infinity}] == 
  Sum[ (-1)^k*Exp[-k], {k, 0, Infinity}] * 
     Integrate[ Exp[-x]/Sin[Pi*x], {x, eps, 1-eps}]

3行目では、整数kに対してSin [Pi *(k + x)] ==(-1)^ k * Sin [Pi*x]が使用されました。

于 2011-04-07T05:34:14.500 に答える
4

私はサーシャに同意する必要があります、積分は収束しているようには見えません。ただし、積分を除外x == 0して分割すると、

Integrate[Exp[-x]/Sin[Pi x], {x, n + 1/2, n + 3/2}, PrincipalValue -> True]

ここn >= 0 && Element[n, Integers]で、交代級数を取得する可能性があるようです

I Sum[ (-1/E)^n, {n, 1, Infinity}] == - I / (1 + E )

今、私はそれをに取り出しただけですがn == 4、それは合理的に見えます。Assumptions -> Element[n, Integers] && n >= 0ただし、上記のMathematicaとの積分では、

If[ 2 n >= 1, - I / E, Integrate[ ... ] ]

これは個々のケースに適合していません。補足として、極が積分領域の境界にある場合、つまり限界が{x, n, n + 1}である場合、sのみが得られますDirectedInfinity。プロットをざっと見てみると、限界がある{x, n, n + 1}場合は厳密に正または負の被積分関数しかないため、無限の値は、補償が不足していることが原因である可能性があり{x, n + 1/2, n + 3/2}ます。で確認し{x, n, n + 2}ますが、評価されていない積分を吐き出すだけです。

于 2011-04-07T16:29:39.997 に答える
4

サイモン、私はあなたの積分にあまり時間を費やしていませんが、固定相近似を調べてみてください。あなたが持っているのは、滑らかな関数(exp)と非常に振動的な関数(sine)です。1/sin(x)関連する作業は、現在、フォームに眉をひそめているところですexp(if(x))

または、次の直列展開を使用することもできますcosecant(極では有効ではありません)。

In[1]:=Series[Csc[x], {x, 0, 5}]
(formatted) Out[1]=1/x + x/6 + 7/360 x^3 + 31/15120 x^5 +O[x]^6

m>-1すべての について、以下があることに注意してください。

In[2]:=Integrate[x^m Exp[-x], {x, 0, Infinity}, Assumptions -> m > -1]
Out[2]=Gamma[1+m]

ただし、ケースを含まない cosecant (ウィキペディアから) の係数を使用してシリーズを合計すると、 に1/x Exp[-x]収束しません[0,Infinity]

c[m_] := (-1)^(m + 1) 2 (2^(2 m - 1) - 1) BernoulliB[2 m]/Factorial[2 m];
Sum[c[m] Gamma[1 + 2 m - 1], {m, 1, Infinity}]

も収束しません...

したがって、積分を無限大に近似できるかどうかはわかりませんが、大きな N までの解に満足している場合は、これらが役立つことを願っています。

于 2011-04-07T08:08:11.587 に答える