0

q,r2 つの変数 ( ) に依存し、余分な積分が 1 つある関数で二重定積分を実行する際に問題があります。ガウス関数で重み付けしたい関数は次のとおりです。

F(q,r)=f(q,r)+int_{0,r}(h(q,r')dr')

そして、ガウスで重み付けするために再び統合する必要があります。

I(q)=int_{0,inf}(F(q,r)^2*g(r)dr)

ガウスg(r)は座標の中心にありますR

ご覧のとおり、主な問題は、配列とスカラーを混在させていることです。ガウスに使用されるのと同じ方法(np.ogridおよび軸の合計)を使用することは解決策になる可能性がありますが、それを実装する方法がわかりません。

import numpy as np
from scipy.integrate import quad
import math as m

R=53.
R0=40.
delta=50.
c=2.
qm, rm = np.ogrid[0.0005:2.0:0.0005, 20:100:500j]

#normalized gauss function
#g(r)
def gauss_grid(r,Rmin,pd):
    def gauss(r,Rmin,pd):
        sigma=1.5
        return (1/sigma)*np.exp(-((r-Rmin)**2)/(2*sigma**2))
    gauss_grid = gauss(r,Rmin,pd)
    #normalization of gaussian
    gauss_grid /= np.sum(gauss_grid)
    return gauss_grid

#spherical function 
#f(q,r) 
def form(q,R):
    return (4/3)*m.pi*3*(np.sin(q*R)-q*R*np.cos(q*R))/(q**3)

#FINAL function
#I(q)
def helfand():
    def F(q,R):
        #integral (0,R) of h(q,r)
        def integral(q,Rmax):
            #h(q,r)
            def integrand(r,q):
                return np.sin(q*r)*(r**2)/(q*r*(1+np.exp(c*(R0-r))))
            return quad(integrand, 0, Rmax, args=(q))[0]
        return (form(q,R)+delta*integral(q,R))**2

    FF_hel=F(qm,rm)
    FF_hel *= gauss_grid(rm,R,pd)
    I=FF_hel.sum(axis=1)
    return I,qm.ravel()

helfand()

*更新* * **

scipy.integrate ライブラリ (を使用quad) を試しましたが、実行できません。q正しい引数 ( ) を次の関数に渡さないようなものです。ここに私がしようとしているものの非常に単純化されたバージョンがあります:

import numpy as np
from scipy.integrate import quad
import matplotlib.pyplot as plt

R=53.
R0=41.
pd=15.
sigma=1.5

def I(q):
    #(function with integral inside) squared
    def FF(q,r):
        def integral_f(q,r):
            def f(r1,q):
                return np.sin(q*r1)
            return quad(f,0,r,args=(q))[0]

        def h(q,r):
            return (r*np.cos(q*r))
        return (h(q,r)+integral_f(q,r))**2

    #gaussian function normalized
    def g(r,R0):
        def gauss(r,R0):
            return (1/sigma)*np.exp(-((r-R0)**2)/(2*sigma**2))
        return gauss(r,R0)/(quad(gauss,0,np.inf,args=(R0))[0])

    #main function to be integrated with gaussian
    def function(r,q):
        return FF(q,r)*g(r,R)

    return quad(function,0,np.inf,args=(q))[0]

q=np.arange(0.001,1.,0.001)
plt.plot(q,I(q))

エラーは言う:

指定された関数は有効な float を返しません。

4

2 に答える 2

0

積分点の境界にまたがる点の単純な 2D 長方形メッシュを作成します。次に、積分を評価するために、各要素よりもガウス求積法を好みます。これは、重み付けされているかどうかにかかわらず、各積分点で関数を呼び出し、直交重みを掛けて、合計することを意味します。

これは、2D 四辺形有限要素と、数値積分による剛性マトリックスの評価に似ています。

SciPy には 2D 求積法があります。自分で書く前にそれを使用します。

http://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.quadrature.html

于 2013-09-27T18:54:41.187 に答える
0

これは 2 つの単積分で計算できると思います。

二重積分を書き出すと、次の 2 つの部分が得られます。

int_{0,inf}(f(q,r)*g(r)dr)+ int_(0,inf)( int_{0,r}(h(q,r')dr')*g(r)博士)

2 番目の積分の順序を交換して取得できます

int_(0,inf)( int_{r',inf}(g(r)dr) * h(q,r')dr')

内部積分は、相補誤差関数で表すことができます。

于 2013-10-02T10:05:35.403 に答える