2

で行われていることに類似したアプローチで、を使用してベクトル値関数の 1 次元積分についてここで与えられたアイデアodeint2dblquad次元に拡張しようとしています。

以下に、私の現在の試みを示します。

import numpy as np
from scipy.integrate import odeint

def _infunc(x, func, gfun, hfun, more_args):
    a = gfun(x)
    b = hfun(x)
    y0 = f(x, a)
    return odeint(lambda v, y: f(x, y, *more_args), y0=y0, t=[a, b] )[1]

def dblodeint(f, a, b, gfun, hfun, args=()):
    y0 = f(a, gfun(a), *args)
    return odeint(lambda v, y: _infunc(y, f, gfun, hfun, args),
                  y0=y0, t=[a, b])[1]

if __name__ == '__main__':    
    def f(x, y):
        return np.array([x*y**2, x**2*y, x**4*y, x**6*y], float)

    def exact_int(a, b, ya, yb):
        return np.array([(b**2-a**2)*(yb**3-ya**3)/6.,
                         (b**3-a**3)*(yb**2-ya**2)/6.,
                         (b**5-a**5)*(yb**2-ya**2)/10.,
                         (b**7-a**7)*(yb**2-ya**2)/14.], float)

    print 'approx:', dblodeint(f, 0, 10, lambda x:0, lambda x:10)

    print 'exact:', exact_int(0, 10, 0, 10)

残念ながら、これは機能していません...次の結果が間違っています:

approx:Repeated convergence failures (perhaps bad Jacobian or tolerances).
Run with full_output = 1 to get quantitative information.
 [ 0.  0.  0.  0.]
exact: [  1.66666667e+04   1.66666667e+04   1.00000000e+06   7.14285714e+07]
4

0 に答える 0