1

出力が数値ではなく新しい関数である Python を使用して関数を統合したいと考えています。たとえば、次の方程式があります (Arnett 1982 から -- 超新星の分析的説明)。

def A(z,tm,tni):
     y=tm/(2*tni)
     tm=8.8             # diffusion parameter
     tni=8.77           # efolding time of Ni56
     return 2*z*np.exp((-2*z*y)+(z**2))

次に、A の積分を見つけて、結果をプロットします。まず、単純に scipy.quad を試しました。

def Arnett(t,z,tm,tni,tco,Mni,Eni,Eco): 
     x=t/tm
     Eni=3.90e+10       # Heating from Ni56 decay
     Eco=6.78e+09       # Heating from Co56 decay
     tni=8.77           # efolding time of Ni56
     tco=111.3          # efolding time of Co56
     tm=8.8             # diffusion parameter 
     f=integrate.quad(A(z,tm,tni),0,x)      #integral of A
     h=integrate.quad(B(z,tm,tni,tco),0,x)  #integral of B
     g=np.exp((-(x/tm)**2))
     return Mni*g*((Eni-Eco)*f+Eco*h)

ここで、B も定義済みの関数です (ここでは示されていません)。A と B はどちらも z の関数ですが、最終的な式は時間 t の関数です。（ここでコードが失敗する原因になっていると思います。）

A と B の積分はゼロから x まで実行されます。x は時間 t の関数です。このままコードを実行しようとすると、「ValueError: 複数の要素を持つ配列の真の値があいまいです。a.any() または a.all() を使用してください」というエラーが表示されます。

そのため、少し検索した後、sympy が適しているのではないかと思いました。しかし、私もこれで失敗しています。

このタスクを完了する方法を教えてください。

どうもありがとう、ザック

4

3 に答える 3

2

A を解析的に積分できます。遅すぎるために何かばかげたことを見逃していないと仮定すると、次のことは役に立ちますか?

import sympy as sy
sys.displayhook = sy.pprint
A, y, z, tm, t, tni = sy.symbols('A, y, z, tm, t, tni')
A = 2*z*sy.exp(-2*z*y + z**2)
expr = sy.integrate(A, (z,0,t)) # patience - this takes a while
expr
# check:
(sy.diff(expr,t).simplify() - A.replace(z,t)).simplify()
# thus, the result:
expr.replace(y,tm/(2*tni)).replace(t,t/tm)

最後の行は、解析形式で A 関数の積分を生成しますが、虚数誤差関数 (scipy.special.erfi() で実行できます) を評価する必要があります。

于 2015-08-06T12:05:23.500 に答える
0

あなたが探しているのはラムダ式だと思います（あなたの言ったことを正しく理解していれば..ラムダ関数の追加情報といくつかの例については、こちらを参照してください）。

彼らができることは、Aで無名関数を定義し、それを返してB関数を取得することです。次のように機能する必要があります。

 def A(parameters):
     return lambda x: x * parameters # for simplicity i applied a multiplication
                                     # but you can apply anything you want to x
 B = A(args)
 x = B(2)

納得のいく回答ができれば幸いです！

于 2015-08-06T11:24:31.553 に答える
0

あなたが得るエラーは、scipy.integrate.quad への間違った呼び出しから来ていると思います:

  • 最初の引数は関数名だけである必要があり、積分はこの関数の最初の変数に対して実行されます。他の変数の値は、args キーワードを介して関数に渡すことができます。
  • scipy.integrate.quad の出力には、積分の値だけでなく、推定誤差も含まれています。したがって、2 つの値のタプルが返されます。

最終的に、次の関数が機能するはずです。

def Arnett(t, z, Mni, tm=8.8, tni=8.77, tco=111.3, Eni=3.90e+10,
           Eco=6.78e+09): 
  x=t/tm
  f,err=integrate.quad(A,0,x,args=(tm,tni))      #integral of A
  h,err=integrate.quad(B,0,x,args=(tm,tni,tco))  #integral of B
  g=np.exp((-(x/tm)**2))
  return Mni*g*((Eni-Eco)*f+Eco*h)

しかし、さらに優れた解決策は、おそらく A と B を分析的に統合してから、ミュリソンが示唆したように式を評価することです。

于 2015-08-06T17:50:35.630 に答える