pythonでLebedev quadratureを使う（更新）

以前にquadpyについて紹介した。
koideforest.hatenadiary.com
しかし、記事が古くなっていて、仕様がいろいろ変更されたようなので、更新する。
詳しい内容は公式のGitHubページを参照。
GitHub - nschloe/quadpy: Numerical integration (quadrature, cubature) in Python

テストコード。

import numpy as np
import quadpy

def f_cart( xyz ):
    # np.shape( xyz ) -> ( 3, number of points )
    n_points = np.shape( xyz )[1]
    return np.ones( n_points )
    # DO NOT WRITE "return 1"

def f_sphere( theta_phi ):
    # np.shape( theta_phi ) -> ( 2, number of points )
    # theta: polar : theta_phi[0]
    # phi: azimuth : theta_phi[1]
    n_points = np.size( theta_phi[0] )
    return np.ones( n_points )
    # DO NOT WRITE "return 1"

#scheme = quadpy.u3._lebedev.lebedev_019()
scheme = quadpy.u3.get_good_scheme(19)

# \int^\pi_0 d\theta \sin\theta \int^{2\pi}_0 d\phi = 4 \pi

xyz_center = ( 0., 0., 0. )
radius = 1.
val_cart = scheme.integrate( f_cart , xyz_center, radius )

val_sphere  = scheme.integrate_spherical( f_sphere )

print( val_cart, val_sphere, 4*np.pi )